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FOREWORD 


This  monograph  provides  a  review  of  the  state  of  the  art  of  modeling  the  important 
processes  of  the  transport  of  mass,  momentum,  and  energy  through  the  lower  atmosphere 
from  the  personal  perspective  of  one  involved  in  turbulent  transport  modeling.  No 
historical  perspective  of  the  vast  literature  published  on  this  subject  is  attempted.  Much  of 
our  understanding  of  atmospheric  boundary-layer  processes  has  come  from  careful 
observation.  For  a  review  of  observational  methods  and  programs,  the  reader  is  referred  to 
the  recent  monograph  by  McBean,  et  al.  (1979).  The  monograph  also  provides  a  very  useful 
review  of  atmospheric-boundary-layer  models.  I  use  this  fact,  plus  the  existence  of  the 
books  by  Brown  (1974),  Plate  (1971),  Wipperman  (1973),  Longhetto  (1980),  and  the  AMS 
workshop  proceedings  edited  by  Haugen  (1973)  and  by  Wyngaard  (1980),  as  justification 
for  further  specialization  of  this  review  to  deal  specifically  with  models  having  a  basic 
theoretical  foundation  in  turbulent  transport  modeling.  Data  and  the  results  of  other  types 
of  models  will  only  be  included  herein  to  the  extent  that  1  feel  they  help  to  place  the  results 
of  turbulent  transport  modeling  in  perspective. 

The  goal  of  this  monograph  is  to  describe  the  approach  of  turbulent  transport 
modeling  to  low-level  atmospheric  dynamics  in  sufficient  detail  to  permit  a  researcher  to 
follow  this  route,  to  review  the  status  of  our  current  understanding  of  boundary-layer 
meteorology,  and  to  describe  some  of  the  applications  where  turbulent  transport  modeling 
can  be  expected  to  have  an  important  impact  in  the  near  future. 

I  wish  to  thank  AGARD  for  this  opportunity  to  organize  my  thoughts  on  modeling 
low-level  atmospheric  turbulence,  and  to  acknowledge  the  following  US  Government 
Agencies;  the  Naval  Air  Systems  Command,  the  Office  of  Naval  Research,  the  Nuclear 
Regulatory  Commission,  the  Environmental  Protection  Agency,  and  the  National 
Aeronautics  and  Space  Administration;  and  the  Electric  Power  Research  Institute  for  their 
original  support  of  my  personal  work  and  that  of  my  A.R.A.P.  colleagues  in  this  area  over 
the  past  few  years.  Naturally,  most  of  the  examples  covered  in  detail  in  this  review  are 
based  on  this  work.  I  apologize  to  other  boundary-layer  modelers  for  not  having  worked 
harder  to  cover  their  work  to  the  same  depth. 
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ABSTRACT 

Y  The  interdependence  between  the  turbulent  transport  of  mass, 
momentum,  and  energy  through  the  lower  levels  of  the  atmosphere,  and  the 
distribution  of  wind,  temperature,  and  species  within  this  layer  are 
reviewed.  Particular  emphasis  is  placed  on  models,  either  analytical  or 
numerical,  which  have  a  basic  theoretical  foundation  in  turbulent 
transport  modeling.  The  results  of  example  model  solutions  are  used  to 
discuss  such  micrometeorological  problems  as:  wind  shear  and  turbulence 
around  airports;  the  prediction  of  low-level  clouds  and  fog;  the 
dispersion  of  industrial  pollutants;  the  Interaction  of  electromagnetic 
radiation  with  turbulent  fluctuations  in  humidity  and  temperature;  and 
the  evolution  of  organized  features  within  the  boundary  layer. 
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Model  coefficient,  see  Eqs.  (2.3.12  and  2.3.13) 

Model  coefficient,  see  Eq.  (2.3.15) 

Mean  species  concentration 

Species  fluctuation 

Specific  heat  at  constant  pressure 

Proportionality  factor  in  Eq.  6.1.10 

General  force  vector  per  unit  mass 

Upward  and  downward  thermal  radiation  fluxes 

Coriolis  parameter  due  to  earth's  rotation 

Gravitational  acceleration 

Mean  value  of  mixing  ratio  of  water  to  air 

Fluctuation  of  water  mixing  ratio 

Mean  and  fluctuating  value  of  liquid  water  mixing  ratio 
Mean  and  fluctuating  value  of  saturated  water  mixing  ratio 
Mean  and  fluctuating  value  of  water  vapor  mixing  ratio 
Measure  of  the  surface  moisture  flux,  wh/U* 

Height  of  vegetation  canopy 

Thermal  conductivity 

Species  diffusivity 

Effective  eddy  thermal  diffusivity 

Effective  eddy  diffusivity  for  moisture 

Coefficients  in  the  moist  air  model,  see  Eqs.  (4.4.6  and  4.4.7) 
Mon in-Obukhov  length 

Eddy  viscosity  length  scale,  see  Eq.  (2.3.1) 

Latent  heat  of  vaporization 
Baroclinic  parameter,  see  Eq.  (4.2.4) 

Mean  pressure 

Pressure  fluctuation 

Heat  source  per  unit  mass 

Species  source  per  unit  mass 

Root-mean  square  velocity  fluctuation 

Universal  gas  constant 

Richardson  nunber 

Gas  constant  for  water  vapor 

Rossby  number 

Radius 

Radius  at  which  maximum  swirl  occurs 

Model  coefficient  In  the  dissipation  of  temperature  variance 

Referenced  background  absolute  teoperature 

Time 

Mean  velocity  components 
Fluctuating  velocity  components 
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W  Mean  vertical  velocity 
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ar  Radiation  diffusivity  factor 
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*<  Earth's  rotation 

“  Frequency  of  pressure  changes  imposed  on  the  boundary  layer 


% 


vi 


1 


i.  Introduction 


The  processes  of  the  transport  of  mass,  momentum  and  energy  through  the  lower  atmosphere  are  of 
extreme  importance  to  a  large  nunber  of  man's  activities.  Perhaps  the  most  basic  necessity  of  life  for 
any  species  is  the  requirement  of  providing  some  means  of  energy  flow  through  the  individual  members  of 
that  species.  A  majority  of  the  species  of  life  on  earth  depend  upon  motions  in  the  lower  levels  of  the 
atmosphere  to  provide  for  at  least  some  part  of  either  the  input  to  the  individual  or  the  disposal  of 
waste  from  it.  This  observation  remains  true  for  a  majority  of  the  machines  invented  by  man.  Thus,  some 
model  of  motions  in  this  part  of  the  atmosphere  is  of  interest  to  a  variety  of  disciplines. 

As  a  result  of  the  widespread  interest,  a  vast  amount  of  literature  has  been  published  on  this 
subject.  No  attempt  will  be  made  to  review,  even  a  reasonable  sampling  of  this  literature.  Rather,  we 
start  by  recognizing  that  it  is  the  small  scale  turbulent  motions  which  generally  play  the  dominant  role 
in  interchanging  mass,  momentum,  or  energy  between  a  surface  and  the  air  above  or  around  it.  The  central 
problem  has  been  to  develop  a  viable  means  of  modeling  turbulent  transport.  Due  to  the  wide  range  of 
scales  involved,  there  is  no  rigorous,  theoretical  solution  for  turbulent  flow.  Turbulent  transport 
models  are  a  blend  of  scaling  based  on  the  physical  conservation  laws  and  empiricism.  Also,  it  appears 
that  the  modeler  frequently  needs  to  apply  a  little  art  to  cover  the  holes  in  the  science. 

The  governing  conservation  equations  for  the  turbulent  transport  of  mass,  momentum,  and  energy  will 
be  discussed  and  the  turbulence  closure  problem  reviewed.  Emphasis  is  placed  on  second-order  closure 
modeling,  because  I  believe  that  this  provides  a  firm  physical  framework  for  extending  the  empirical 
information  gained  from  select  experiments  to  much  more  general  problems.  While  it  is  true  that  much, 
perhaps  even  most,  of  our  present  understanding  of  dynamics  in  the  atmospheric  boundary  layer  has  come 
from  simpler  scaling  laws  applied  to  reliable  field  data,  I  believe  that  second-order  closure  has  in 
recent  years  demonstrated  that  it  will  have  a  major  impact  on  solving  those  problems  where  the  turbulent 
fluxes  or  turbulent  variances  are  the  quantities  of  most  interest.  A  rough  review  of  the  scales  of  the 
boundary-layer  thickness  over  which  the  strongest  interaction  between  the  surface  and  the  atmosphere  occur 
is  also  given  in  Chapter  2.  This  monograph  is  entitled  "Modeling  the  Lowest  1  Km  of  the  Atmosphere",  to 
emphasize  that  this  is  the  region  of  most  concern  to  us.  However,  the  scaling  analysis  in  the  next 
chapter  makes  it  evident  that  the  turbulent  boundary  layer  generated  from  the  surface  may  exceed  1  Km  in 
thickness  under  certain  conditions  and  be  much  less  than  1  Km  under  other  conditions.  In  fact,  under  the 
conditions  most  often  occurring  in  nature,  there  may  be  a  number  of  layers  governed  by  different  dynamical 
behaviors  superimposed  at  any  one  time.  In  discussing  these  idealized  layers,  we  will  deal  with  the 
complete  layer,  but  keep  in  mind  that  we  are  most  interested  in  that  portion  of  the  layer  which  falls 
within  1  Km  height. 

Next,  we  look  in  detail  at  some  of  the  separate  balances  which  can  permit  ideal  similarity  laws  to 
exist  in  at  least  a  portion  of  the  lowest  1  Km  altitude.  The  first  of  these  is  the  surface  layer  which 
may  only  extend  up  a  few  meters  or  as  much  as  a  few  hundred  meters.  It  is  these  surface-layer, 
Monin-Obukhov,  similarity  laws  which  provide  the  most  reliable  information  in  current  widespread  use  in 
micrometeorology.  The  assumptions  necessary  for  validity  of  the  surface  layer  functions  are  often 
fulfilled.  This  is  not  so  often  true  of  the  other  ideally  similar  layers.  However,  they  do  form  a  basis 
for  our  understanding  the  different  types  of  dynamics  which  the  lower  atmosphere  will  support.  In 
particular,  they  clearly  show  the  role  of  the  coriolis  force  in  turning  the  wind  vector  through  a 
significant  angle  within  this  height  domain,  the  role  of  heat  transfer  between  the  surface  and  the 
atmosphere  in  either  amplifying  or  damping  the  turbulence  generated  by  the  moving  wind  over  the  stationary 
surface,  and  the  role  of  the  structure  of  the  surface  canopy  In  defining  the  detailed  interaction  between 
the  surface  and  the  atmosphere.  From  these  ideal  similarity  laws  it  is  clear  that,  not  only  do  the  wind, 
temperature,  and  species  control  the  interchange  of  mass,  momentum,  and  energy  between  the  surface  and 
atmosphere,  but  it  is  just  as  valid  to  view  the  surface  fluxes  of  mass  (in  the  form  of  water  vapor), 
momentum,  and  energy  as  predominantly  controlling  the  motions  of  the  atmosphere  in  the  lowest  1  Km.  These 
similarity  laws  have  formed  the  basis  for  most  of  the  currently  used  engineering  methods  for  estimating 
such  things  as:  the  surface  drag,  the  angle  the  wind  at  the  surface  departs  from  the  direction  of  the 
wind  at  the  top  of  the  boundary  layer,  wind  shear  within  this  layer,  surface  heat  transfer,  and  water 
evaporation  from  the  surface. 

In  the  real  world,  the  interaction  of  turbulent  transfer  with  other  physical  mechanisms  modifies  the 
limiting  similarity  solutions.  Typical  time  variations  of  the  surface  heat  transfer  caused  by  the  diurnal 
solar  radiation  delivered  to  the  land  drive  the  boundary  layer  between  the  strongly  stable  nocturnal  layer 
and  the  buoyantly  driven,  well-mixed  layer  of  the  afternoon.  A  strong  hysteresis  effect  exists  with  the 
morning  “neutral"  boundary  layer  quite  distinct  from  the  evening  “neutral"  layer.  Other  effects,  such  as 
small  horizontal  variations  in  temperature,  of  the  order  of  1!C  per  100  Km,  lead  to  important  baroclinic 
gradients  in  the  wind.  Example  calculations  showing  the  Influence  of  typical  unsteady  and  spatial 
variability  demonstrate  that  the  steady  state,  spatially  homogeneous,  neutral  planetary  boundary  layer 
seldom  comes  close  to  being  realized  in  nature. 

The  lowest  1  Km  of  the  atmosphere  plays  a  key  role  in  controlling  the  cycling  of  water  through  the 
atmosphere.  On  a  global  average  basis,  the  turbulent  transport  of  water  vapor  through  this  atmospheric 
layer  must  balance  the  average  precipitation  of  liquid  water,  about  1  m  per  year  (Coantic,  197d).  While 
It  is  the  small  scale  turbulent  eddies  which  drive  the  transport  of  water  vapor  through  the  boundary 
layer,  it  is  the  release  of  latent  energy  transported  with  the  water  vapor  which  is  responsible  for 
fueling  most  of  the  larger  scale  atmospheric  motions.  These  larger,  synoptic  scale  meteorological  motions 
are  not  analyzed  herein,  but  instead  are  treated  as  boundary  conditions  at  the  top  of  the  boundary  layer. 
Latent  energy  release  within  the  boundary  layer  in  the  form  of  low-level  clouds  or  fog  also  directly 
affects  the  structure  of  the  wind  and  temperature  in  this  layer. 

An  Indirect,  but  often  even  stronger.  Influence  of  condensed  water  vapor  In  the  boundary  layer  is  its 
interaction  with  thermal  radiation.  On  the  global  scale,  the  energy  to  drive  the  atmosphere  ultimately 
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traces  back  to  the  sun,  with  thermal  radiation  the  principal  mode  of  energy  transfer  through  tne 
atmosphere.  Under  clear  conditions,  the  relative  transparency  of  the  atmosphere  allows  this  radiant 
energy  to  pass  through  the  lowest  1  Km  without  significant  interaction.  However,  the  presence  of 
condensed  water  droplets  in  the  air  greatly  reduces  thermal  radiation  absorption  depths  and  permits  a 
strong  interaction  between  turbulent  transport  and  thermal  radiation.  Thus,  the  influence  of  thermal 
radiation  must  be  incorporated  in  any  boundary-layer  model  which  seeks  to  simulate  low-level  cloudiness  or 
fog. 


The  vertical  variations  normal  to  the  surface  are  generally  of  most  concern  in  modeling  fluxes 
between  the  surface  and  the  atmosphere,  but  the  boundary  layer  contains  a  rich  variety  of 
three-dimensional  flow  features  which  are  of  increasing  Interest.  In  simulating  the  boundary  layer  it  is 
necessary  to  decide  when  to  seek  to  resolve  particular  features  rather  than  to  ensemble  average  over  them 
as  part  of  the  turbulence.  Chapter  S  considers  a  number  of  particular  features  involving  two-dimensional 
symmetry,  either  planar  or  axisymmetric.  Simulations  of  these  features  are  computed  by  incorporating  the 
boundary  conditions  for  this  particular  2-0  feature  into  the  boundary  conditions  for  the  mean  variable 
while  representing  the  smaller  scale  three-dimensional  eddies  with  the  turbulence  closure  modeling.  The 
features  simulated  vary  from  a  simple  step  change  in  surface  roughness  conditions,  to  the  most  violent 
eddy  which  occurs  in  nature,  the  tornado. 

The  simulation  of  the  common  longitudinal  roll  vortices  which  extend  the  full  height  of  the  boundary 
layer  provide  an  opportunity  to  examine  the  relative  contributions  of  the  smaller  scale  random  turbulence 
and  the  organized  large  scale  rolls.  Both  play  an  important  role  in  transporting  mass,  momentum,  and 
energy  through  the  lowest  1  Km  of  the  atmosphere. 

In  particular  engineering  applications,  it  is  generally  desirable  to  work  with  the  simplest  model 
which  provides  a  reliable  answer  to  the  question  at  hand.  In  Chapter  6,  several  particular  applications 
are  discussed  which  provide  examples  where  the  turbulent  fluxes  and  variances  are  the  quantities  of  most 
concern.  This  provides  a  sampling  of  areas  where  the  relatively  complex  turbulent  transport  modeling 
reviewed  herein  is  currently  making  contributions  to  engineering  problems.  The  list  is  not  intended  to  be 
all  inclusive,  rather  an  attempt  was  made  to  describe  each  individual  application  in  sufficient  detail  to 
make  it  useful  to  the  interested  reader. 
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II.  Model  Formulation 


2.1.  Basic  Equations 

The  basic  equations  governing  the  transport  of  mass,  momentum,  and  energy  in  the  atmospheric  boundary 
layer  have  been  well  known  for  some  time.  They  may  be  written  as: 


Mass  continuity: 


Momentum: 


Energy: 


Species: 
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The  Incompressible  form  of  the  continuity  equation  is  appropriate  for  the  atmosphere  as  long  as  we 
restrict  attention  to  velocities  much  less  than  the  speed  of  sound  and  to  shallow  convection,  l.e., 
vertical  motions  which  extend  for  a  height  ouch  less  than  the  scale  height  of  the  atmosphere.  The  most 
Important  effects  of  density  changes  In  the  boundary  layer  are  those  resulting  In  combination  with  gravity 
acceleration  as  first  pointed  out  by  Boussinesq  (1903).  This  buoyant  force  term  and  the  coriolis  body 
force  due  to  the  earth's  rotation  are  the  primary  terms  which  will  be  used  in  F^.  The  source  terms  (JL  and 
Qr  in  the  energy  and  species  equations  will  generally  be  zero  except  for  that  due  to  thermal  radiation 
flux  divergence  in  Q9. 

The  temperature  variable  Is  chosen  as  the  potential  temperature 


e  »  T  -  Tref  ♦  z 
CP 


(2.1.5) 


so  that  the  temperature  changes  experienced  by  adiabatic  expansion  or  compression  accompanying  vertical 
motions  will  be  automatically  Incorporated. 

The  primary  difficulty  with  Eqs.  (2.1.1  -  2.1.4)  Is  that  although  these  equations  may  be  solved 
rather  straightforwardly  with  today's  computing  resources  for  simple  boundary  conditions,  the  large  range 
of  scales  Involved  In  the  motions  of  the  atmospheric  boundary  layer  prohibits  a  complete  detailed 
solution.  Atmospheric  motions  range  from  small  dissipative  eddies,  as  small  as  10~3  meters,  to  the 
synoptic  weather  features  of  order  10*  meters.  It  appears  essential  to  average  over,  at  least,  the  lower 
end  of  this  spectrum  of  scales.  These  smallest  scale  motions  tend  to  be  random  in  character  and  thus  most 
susceptible  to  statistical  analyses. 


When  Eqs.  (2.1.1  -  2.1.4)  are  averaged  over  an  ensemble  of  realizations  the  equations  may  be  written 
as: 


Mass  continuity: 
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The  dependent  variables  In  Eqs.  (2.1.6  -  2.1.9)  have  been  split  into  an  ensemble  mean  value,  denoted  by 
capitalization,  representing  the  scales  of  motion  to  be  resolved,  and  a  fluctuating  value,  denoted  by 
lower  case  symbols,  representing  the  unresolved  motions. 

The  price  for  averaging  Eqs.  (2.1.1  -  2.1.4)  to  reduce  the  range  of  scales  to  be  resolved  is  high. 
Equations  (2.1.6  •  2.1.9)  do  not  form  a  complete  set  mtil  the  second-order  correlations,  u,uj,  u,«,  and 
u,c  and  u,c,  can  be  determined  or  related  to  the  mean  variables.  This  is  the  central  closure  problem  of 
turbulent  transport  theory.  The  larger  the  scales  of  motion  included  in  the  unresolved  averages,  the  more 
sensitive  the  task  of  correctly  relating  the  second-order  correlations  to  the  resolved  flow  field.  On  the 
other  hand,  the  smaller  the  scales  of  motion  included  within  the  resolved  motion,  the  more  difficult  the 
task  of  solving  for  this  mean  motion.  The  choice  of  the  break  between  the  resolved  and  unresolved  scales 
calls  for  a  compromise  dictated  by  the  particular  application. 

A  closure  scheme  with  general  validity  is  highly  desirable.  A  method  of  incorporating  additional 
physical  conservation  relationships  into  the  closure  scheme  is  to  begin  with  the  equations  for  u,uj,  u,u, 
and  u,c  derived  from  Eqs.  (2.1.1  -  2.1.4).  With  the  aid  of  Eqs.  (2.1.6  -  2.1.9)  and  a  fair  amount  of 
algebraic  manipulation,  these  may  be  written  as: 
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(2.1.12) 


As  proposed  by  Lilly  (1969)  and  Oonaldson  (1972,  1973),  the  physical  constraints  Inherent  within 
Eqs.  (2.1.10  -  2.1.12)  can  provide  strong  guidelines  for  modeling  turbulent  transport  in  the  atmosphere 
assuming  the  terms  involving  higher-order  correlations  can  be  represented  in  a  consistent  manner.  A 
variety  of  models  have  been  presented  by  researchers  for  closing  the  equations  at  this  second  order. 
Before  proceeding  to  some  of  these  models,  let  us  discuss  the  scaling  inherent  in  Eqs.  (2.1.10  -  2.1.12) 
in  order  to  outline  the  dependency  of  the  boundary-layer  thickness  on  the  pertinent  variables. 


2.2.  Approximate  Scaling  Relationships 

When  characteristic  lengths  L  and  S,  and  characteristic  velocities  V  and  W  are  introduced  for  the 
horizontal  and  vertical  components  respectively,  mass  continuity,  Eq.  (2.1.6)  shows  that 


W 


(2.2.1) 


The  further  introduction  of  F  for  the  characteristic  body  force,  q  for  the  characteristic  turbulence 
velocity,  and  t  for  the  characteristic  time,  permits  the  momentum  equation,  Eq.  (2.1.7),  to  be  scaled  as 


V  +  V2 
t  L 


(2.2.2) 


To  help  Interpret  this  equation,  it  is  advantageous  to  add  the  Reynolds  stress  equation,  Eq.  (2.1.10) 
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A  characteristic  length  a  for  the  turbulence  has  been  introduced  into  Eq.  (2-2.3) -  It  is  reasonable  to 
expect  a  to  be  some  fraction  of  6. 


Equations  (2.2.2)  and  (2.2.3)  lead  to  a  number  of  different  possibilities.  First,  if  the  flow  is 
steady  and  horizontally  homogeneous,  i.e.,  both  z  and  L  ♦  «,  then  from  Eq.  (2.2.2) 


& 


(2.2.4) 


and  from  Eq.  (2.2.3) 
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(2.2.5) 


If  coriolis  body  forces  are  considered  then  F  ~  fV,  and  thus 


6 


(2.2.6) 


For  typical  values  of  (a/6)  ~  0.1,  V  ~  10  m/sec,  f  ~  10-4sec-1,  this  leads  to  6  ~  1  Km  for  the  scale  of 
the  planetary  boundary  layer  (PBL)  under  steady  state,  constant  density,  spatially  homogeneous  conditions. 
These  ideal  conditions  seldom  apply  in  the  atmosphere. 

Unsteady  variations  in  F  may  be  expected  to  reduce  6.  For  example,  if  the  geostrophic  wind  is 
assumed  to  have  a  periodic  variation  with  frequency  u,  then  the  time  scale  should  scale  with  w'1  and 
Eq.  (2.2.6)  may  be  modified  to 
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(2.2.7) 


Thus  the  boundary  layer  associated  with  mesoscale  features  with  time  scales  less  than  s:  3  hours,  such  as 
convective  storms,  should  be  expected  to  be  thinner  than  the  normal  steady  state  values. 


Small  differences  in  air  density  also  play  a  very  important  role  through  Fj  or  F. 
be  order  g66/T0  while  F  remains  order  fV,  then  it  appears  possible  to  have 


If  F j  is  taken  to 
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(2.2.8) 


It  will  be  seen  in  Section  3.4  that  this  is  appropriate  for  an  equilibrium  layer  stabilized  by  negative 
buoyancy.  In  practice,  this  limit  also  appears  to  be  only  rarely  reached.  Instead,  it  is  possible  for 
Eq.  (2.2.3)  to  balance  on  a  shorter  time  scale  than  f"i.  In  the  stabilizing  case,  the  boundary  layer 
thickness  determined  by  a  balance  of  the  3  terms  on  the  right-hand  side  of  Eq.  (2.2.3)  is 


A  V2T0 
6  gas 


(2.2.9) 


In  the  unstable  case,  the  last  two  terms  in  Eq.  (2.2.3)  attempt  to  balance,  but  if  the  lapse  rate,  *e/»z, 
remains  negative,  the  turbulence  will  continue  to  grow  unbounded.  Under  these  conditions,  Eq.  (2.2.3)  may 
be  interpreted  as  only  yielding  the  scale  of  q2  for  any  6, 
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with  6  growing  compatibly  with  the  mean  potential  temperature  equation,  Eq.  (2.1.8),  or  more  frequently 
limited  by  a  capping  stable  temperature  lapse  rate. 

When  Eq.  (2.2.9)  is  Interpreted  in  terms  of  lapse  rate,  r  *  3e/»z,  then 
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(2.2.11) 


Relatively  small  stable  lapse  rates  (i,  I0~3°c/m)  lead  Eq.  (2.2.11)  to  yield  values  of  6  much  lower  than 
those  obtained  from  Eq.  (2.2.6).  Under  light  wind  conditions  (~  1  m/sec)  and  strong  lapse  rates 
(~  10*loC/m),  as  may  be  set  up  In  the  early  morning  hours  after  a  clear  night  of  surface  radiational 
cooling,  4  can  be  as  small  as  ~  10  meters. 

Spatial  homogeneity,  although  still  infrequent,  is  more  readily  obtained.  It  requires  that 
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There  are  many  places  on  earth  where  the  surface  Is  relatively  homogeneous  over  distances  of  the  order  of 
100  Km.  A  frequent  obstruction  to  spatial  homogeneity  comes  not  from  the  surface  roughness  boundary 
conditions,  but  from  thermal  boundary  conditions.  Horizontal  temperature  gradients  of  the  order  of  1°C 
per  100  Km,  or  larger,  lead  to  vertical  variations  in  the  pressure  gradient  driving  the  PBL  which  can 
significantly  influence  the  flow.  Likewise,  horizontal  variations  in  the  free  stream  wind  of  the  order  of 
1  m/ sec  per  100  Km  which  lead  to  horizontal  divergence  of  the  order  of  lO'Ssec*1  lead  to  vertical 
velocities  of  the  order  of  1  cm/sec  at  the  top  of  the  boundary  layer.  Although  seemingly  small,  these 
vertical  velocities  are  sufficiently  large  to  have  significant  effects  on  the  vertical  gradients  of 
velocity,  temperature,  and  humidity  at  the  top  of  the  boundary  layer. 

As  a  way  of  sunmarizing  the  boundary  layer  scaling,  it  is  possible  to  combine  the  previous  equations 
in  an  approximate  fashion  and  write 


_ (a/6)2V _ 
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(2.2.13) 


Whichever  of  the  frequencies  is  largest;  the  coriot  is  f.  the  unsteadiness  in  the  geostrophic  winds  ",  the 
frequency  of  the  stable  oscillations,  (A/6)3/2(Yg/T0)172,  related  to  the  Brunt-Vaisala  frequency,  or  the 
frequency  associated  with  the  spatial  inhomogeneity,  V/L,  will  control  the  boundary  layer  thickness.  If 
the  density  stratification  is  locally  unstable,  then  6  will  increase  until  it  is  capped  by  a  stable  layer. 
Although  the  largest  frequency  may  be  expected  to  control  the  thickness  over  which  the  strongest 
adjustments  between  surface  and  free  stream  conditions  occur,  in  nature  we  can  expect  these  thin 
adjustment  layers  to  coexist  in  the  form  of  young  shallow  boundary  layers  occurring  within  older  deeper 
layers. 


2.3.  Turbulence  Closure 

A  complete  review  of  the  various  approximations  used  by  different  researchers  to  close  the  set  of 
equations  given  in  Section  2.1  is  beyond  the  scope  of  this  monograph.  Some  reviews  of  this  subject 
appearing  In  the  literature  are  those  by  Mellor  and  Herring  (1973),  Reynolds  (1976),  Lewellen  (1977), 
Lumley  (1979,  1980),  Launder,  et  al.  (1975),  and  Zeman  (1981).  It  currently  remains  an  area  of  active 
research.  What  I  will  attempt  to  do  is  provide  a  flavor  of  some  of  the  different  approaches  and 
comproml ses  being  pursued  and  provide  some  guide  posts  to  the  literature  for  the  reader  who  wishes  to 
delve  more  deeply  Into  this  subject. 

The  basic  compromise  which  must  be  made  is  between  complexity  to  incorporate  as  true  a  representation 
of  the  physics  Inherent  in  Eqs.  (2.1.10  -  2.1.11)  as  possible  and  simplicity  to  make  model  calculations 
possible  for  more  problems  of  Interest.  It  is,  of  course,  not  possible  to  equate  complexity  with  physical 
rigor,  but  It  does  appear  fairly  universally  true  that  practical  models  require  a  liberal  use  of  judicious 
Simplifying  assumptions.  The  simplest  models  we  wish  to  consider  here  are  those  that,  at  least,  make  some 
use  of  the  physics  involved  in  the  Reynolds  stress  equations  in  parameterizing  the  flux  relationships  used 
to  close  the  mean  equations. 

A  table  covering  the  wide  variety  of  mixing  length  models  that  have  been  applied  to  the  PBL  problem 
Is  given  in  McBean,  et  al.  (1979).  The  most  popular  of  these  involves  formulae  of  the  type 
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(2.3.1) 


where  t(z/4)  and  f tRI  *  C@g39/az]/[(au/3z)2  +  (av/8z)2])  are  analytical  fits  to  semi -empirical 
parameterizations.  The  conditions  necessary  for  such  a  form  to  be  consistent  with  Eq.  (2.1.10)  are 
interesting.  First,  it  requires  assuming  the  turbulence  to  be  in  local  equilibrium  so  that  the  left  hand 
side  of  the  equation  as  well  as  the  diffusion  contributions  of  the  higher-order  terms  on  the  right  hand 
side  are  negligible.  Second,  it  requires  assuming  an  expansion  about  isotropic  turbulent  conditions.  In 
spite  of  the  fact  that  these  two  assumptions  are  somewhat  contradictory,  judicious  choices  of  t  and  f  can 
permit  this  formulation  to  lead  to  useful  results  for  many  problems. 
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At  the  other  end  of  the  spectrum  are  the  models  (e.g.,  Andre,  et  al..  1976;  Lumley,  Zeman,  and 
Sless,  1978)  which  carry  the  closure  to  third  order  to  more  correctly  model  the  third-order  terms 
appearing  In  Eq.  (2.1.10).  Such  an  approach  should  be  expected  to  yield  Improvements  In  the  dynamics  of 
the  turbulence,  particularly  In  the  representation  of  the  third-order  turbulent  correlations.  Our 
approach  here  is  to  describe  a  compromise  which  imposes  a  relatively  simple  closure  approximation  to  allow 
Eqs.  (2.1.10  -  2.1.12)  to  be  carried.  With  this  choice,  we  expect  the  computational  results  for  the 
appropriate  second-order  correlations  to  approximate  the  observations  of  such  quantities  as  the  Reynolds 
stress,  the  heat  flux,  species  flux  and/or  temperature  variance  for  which  we  seek  an  engineering  answer. 
We  do  not  expect  all  the  third-order  correlations  to  be  very  well  represented.  With  this  compromise  in 
mind,  let  us  look  at  simple  models  for  the  higher  terms  in  Eqs.  (2.1.10  -  2.1.12). 

The  terms  which  must  be  modeled  may  be  grouped  as  dissipation  terms: 
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terms  Involving 


Interaction  with  the  fluctuating  pressure  gradient: 
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and  third-order  correlations: 
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Lumley  (1979)  views  second-order-closure  modeling  as  assuming  that  the  length  and  time  scales  of  the 
turbulence  are  small  relative  to  the  length  and  time  scales  characterizing  the  distributions  of  the  mean 
quantities;  i.e.,  that  the  flowfield  is  quasisteady  and  quasi  homogeneous.  Although  this  assumption 
permits  the  development  of  a  consistent  expansion  procedure,  it  imposes  a  limitation  which  does  not  appear 
to  be  either  necessary  or  desirable.  It  is  exactly  those  Situations  where  the  length  and  time  scales  of 
both  the  turbulence  and  the  mean  flow  are  of  the  same  order  where  first-order  closure  is  most  seriously  in 
error  and  thus  for  which  we  look  most  urgently  to  second-order  modeling  for  help.  Since  the  limitation  is 
only  imposed  if  it  is  introduced  as  part  of  the  technique  for  chosing  the  modeled  terms,  it  seems  to 
provide  strong  motivation  not  to  rely  too  heavily  on  this  otherwise  elegant  expansion  procedure  for 
determining  model  coefficients. 

Dissipation  terms.  The  dissipation  terms  measure  the  effect  of  molecular  smoothing  on  the  structure 
of  the  turbulent  correlations.  Even  in  high-Reynolds-number  flow,  we  expect  viscous  dissipation  to  be  the 
major  loss  mechanism  for  turbulent  kinetic  energy.  Owing  to  the  nonlinear  terms  in  the  Navier-Stokes 
equation,  a  reduction  in  viscosity  is  compensated  by  a  reduction  In  scale  of  the  smallest  eddies  in  the 
flow.  The  dissipation  eddies  are  much  smaller  than  the  eddies  that  receive  their  energy  directly  from  the 
mean  flow  when  the  Reynolds  mmber  is  large,  and  one  may  assune  that  they  are  statistically  independent  of 
the  mean- flow  geometry.  • 

Just  as  the  end  of  the  cascade  contains  no  information  on  the  scale  of  the  large  eddies,  the  breakup 
rate  of  the  large  eddies  should  be  independent  of  v.  Therefore,  for  high  Reynolds  mmber,  it  appears 
dimensionally  correct  to  have 
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where  q  is  the  root -mean-square  value  nr  the  total  velocity  fluctuation  and  a  is  a  macroscale  of  the 
eddies.  This  corresponds  to  taking  the  Taylor  microscale  x  proportional  to  A/(a  +  bqa/v)1'2,  where  a  and 
b  are  constants,  as  suggested  by  Rotta  (1951).  Since  the  viscous  dissipation  process  is  an  Isotropic 
process,  most  Investigators  model  the  dissipation  terms  as  an  isotropic  term,  i.e.. 
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An  additional  anisotropic  term  seems  appropriate  when  low  Reynolds  numbers  are  permitted  (lewellen,  1977), 
but  since  the  Reynolds  number  in  the  atmospheric  boundary  layer  is  almost  always  large,  we  will  neglect  it 


8 


The  corresponding  models  for  £9  and  ec  are 


j2Ul  s2*) 

v6— y  keu  -r-2 


38  38 

3Xk  3xk 


.  *C  3c 
C  3*k  *k 


Cc  A,  > 


where  Ay  and  Ac  need  not  equal  A. 

Several  Investigators  prefer  to  calculate  c  from  a  dynamic  equation  obtained  by  modeling  its 
governing  equation.  This  will  be  discussed  later  in  connection  with  the  scale  equation  since,  through 
Eq.  (2.3.5),  it  is  equivalent  to  an  equation  for  A. 

Pressure  Correlations.  The  correlations  involving  pressure  fluctuations  redistribute  the  turbulent 
energy  produced  by  the  production  terms.  This  can  perhaps  be  seen  more  readily  by  rewriting  the 
velocity-pressure  correlation  as 


ill  JL 

P 


+  £iJ£. 

P  3*1 


.  £  /i£i  +  i£i 

P\3Xj  3Xi 


(2.3.7) 


Since  the  flow  is  incompressible,  there  will  be  no  contribution  from  the  last  term  in  Eq.  (2.3.7)  to  the 
kinetic  energy  of  the  turbulence,  q2/2  =  u^u^/2.  It  contributes  only  to  a  redistribution  of  energy 
between  the  Reynolds  stresses.  A  volume  integral  of  the  first  two  terms  on  the  right-hand  side  of 
Eq.  (2.3.7)  over  any  finite  region  of  turbulence  bounded  by  laminar  flow  will  yield  zero.  Thus  these  two 
terms  tend  to  contribute  to  a  spatial  redistribution,  i.e.,  are  diffusion  terms. 

All  investigators  appear  to  follow  Rotta  (1951)  in  modeling  at  least  one  contribution  of  the  pressure 
correlation  as  a  return-to- isotropy  term  proportional  to  the  extent  to  which  the  flow  is  anisotropic, 
i.e.. 


“(“("j-'u  f) 


(2.3.8) 


The  term  in  parentheses,  the  departure  from  isotropy,  provides  the  correct  tensor  symmetry  and  q/A 
provides  the  correct  dimensionality.  Donaldson  (1972,  1973)  adds  to  this  a  spatial  diffusion  term  with 
proper  symmetry  and  dimensionality  to  represent  the  first  two  terms  on  the  right-hand  side  of  Eq.  (2.3.7). 


C3q  / qA  [  3  /  3ukUj  \  3  /  3ukUj  \ 

t  p  - 7) +  c>  p  r  7?) + ^  r ) 


(2.3.9) 


Most  Investigators  add  terms  simulating  the  pressure  contributions  in  the  interactions  between  mean  shear 
and  turbulence,  and  between  buoyancy  and  turbulence.  As  given  by  Zeman  (1981)  these  are 


’*'» j (b)  °  '  2<12  jsij  +  y  ^Sikbkj  +  Sjkbki  -  |  SbAij^J 


(2.3.10) 


8U1  +  8Uj  -  ~  8uk 


(2.3.11) 


1  /i£i  + 

2  ysxj  3X1 J 


£i£i  1  * 

q2  '  7  4lj 


where 
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The  coefficients  In  Eq.  (2.3.10)  were  determined  by  matching  to  so-called  rapid-distortion  calculations 
Involving  a  linear  expansion  around  an  Isotropic  state.  Unfortunately  this  model  is  inconsistent  with  the 
majority  of  boundary-layer  data  (Zeman,  1981)  and  various  authors  have  altered  the  constants  in  some  way 
to  achieve  agreement  with  the  data  of  most  interest  to  them.  Lumley  (1979)  argues  that  the  basic  problem 
is  that  the  constant  C3,  in  Eq.  (2.3.8)  is  really  a  function  of  by  by.  Another  approach  for  boundary 
layer  flows  Is  to  add  near-surface  corrective  functions  representing  the  contribution  of  the  surface 
effect  on  fluctuating  pressure  terms  (Gibson  and  Launder,  1978). 

It  is  not  surprising  that  these  pressure  terms  show  a  lot  of  variation  for  different  flows,  since 
these  pressure  terms  involve  spatial  integrals  over  the  domain  of  the  pressure  variation.  Their 
parameterization  in  terms  of  local  properties  of  the  mean  flow  variables  and  second  order  turbulent 
correlations  must  allow  for  the  complete  myriad  of  turbulent  structures  possible.  No  “standard" 
treatments  for  these  terms  are  available  at  present.  Several  of  the  calculations  exemplified  In  the  later 
chapters  will  Involve  setting  the  terms  in  Eqs.  (2.3.10  -  2.3.11)  to  zero  following  Lewellen  (1977)  and 
Andre,  et  al.  (1978).  Such  a  model  involves  a  compromise  between  known  deficiencies  In  certain  limits  and 
the  simplicity  of  allowing  the  model  to  remain  Invariant  over  a  wide  range  of  applications.  Also  the 

to  provide  a  model  which  always 
<  uaua  luug,  and  let  |uaup|  >  0] 
ese  conditions. 

If  only  the  leading  return-to- isotropy  term  is  retained  the  corresponding  tenperature  and  species 
pressure  terms  may  be  modeled  as 


wiuuinakiun  ut  tne  simple  Terms,  tqs.  u.u.o  ana  z.j.b;  can  be  snown 
satisfies  the  real izabil ity  relations  [uau«  >  0  for  a  =  p,  uau«2 
(Schumann,  1977),  while  the  more  complex  forms  do  not  always  satisfy  th 


_e JL  hi _ 

P  3X,  *  Ae  “i6 


(2.3.12) 


c  «p 

P  3Xi 


Aq - 

^UiC 


(2.3.13) 


Third-Order  Velocity  Correlations.  These  terms  represent  a  process  by  which  the  turbulent 
correlation  is  transferred  from  one  part  of  the  flow  to  another  without  any  net  production  or  loss.  This 
can  be  demonstrated,  as  it  was  for  the  pressure  diffusion  term  of  Eq.  (2.3.7),  by  integration  over  a 
finite  volume.  When  the  volume  is  bounded  by  either  laminar  or  by  homogeneous  isotropic  flow,  the 
integral  will  vanish.  The  most  popular  modeling  of  this  transport  term  is  as  gradient  diffusion,  although 
Bradshaw  (1972)  has  suggested  that  it  could  well  be  algebraic. 

A  number  of  different  gradient  diffusion  forms  have  been  used  as a_  model  for  this  term. 

Donaldson  (1972)  proposed  a  form  that  satisfies  the  tensor  symmetry  of  u^u,  with  a  scalar  diffusion 
coefficient. 


3  _ 

°1j  *  —  tURUiUj) 

3*k 


Hanjalic  and  Launder  (1972)  used  a  form  with  tensorial 


diffusivity 


+ 


3»ku1  +  3UkUj  \ 
axj  ax-)  J 


They  obtained  this  form  by  a  "firm  pruning"  of  the  exact  equation  for  utUfUj.  Others  have  Included  more 
Information  from  the  third  order  equations  as  reviewed  by  Zeman  (1981).  Andre,  et  al.  (1978)  carries  21 
rate  equations  for  third  order  correlations  in  their  simulation  of  the  atmospheric  boundary  layer. 


The  simplest  possible  form  for  all  of  the  diffusion  terms  in  Eq.  (2.3.4)  is  given  by 


a 

3*j 


>3 


3*j  [qA  3*j 


>] 


(2.3.14) 


Although  this  form  satisfies  the  tensor  symmetry  of  Dy,  it  does  not  satisfy  the  synmetry  of  Uku^uj 
Itself. 

Modeled  Equations.  Our  philosophy  is  to  first  attempt  calculations  with  the  simplest  possible 
second-order  closure  models.  We  choose  to  use  Eq.  (2.3.6)  for  the  dissipation  term,  Eq.  (2.3.8)  for  the 
pressure-velocity  interaction,  and  the  pressure  contribution  to  diffusion  Incorporated  with  the  velocity 
diffusion  for  one  diffusion  term  like  Eq.  (2.3.14): 


V*M  jili 

3*k  3*k  3A 


(2.3.15) 
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(2.3.16) 


(2.3.17) 


In  this  restatement  of  the  modeled  terms,  C3  has  been  set  equal  to  unity  by  using  Eq.  (2.3.8)  to  define 
the  macroscale  A.  The  other  coefficients  have  been  assigned  the  symbols  given  by  lewellen  (1977). 

When  similar  terms  are  used  in  the  heat-flux,  species-flux,  temperature-variance,  and 
species-variance  equations,  the  modeled  set  of  second-order  correlation  equations  may  be  written  as 


at 


_ JUj  _ aUi  _  _ 

-  uiuk  ^  *  ujuk~  +  MlUJB  + 


a 


+  *c 


a*k 


q/ — 

-  «1j 


(2.3.18) 


auiti  auitt 

-  ♦  Uj  — - 

at  J  axj 


-  ae  —  aui 

-  ujUj -  -  Ujd  - —  +  t  gj»*  + 

axj  axj 


(2.3.19) 


(2.3.20) 


(2.3.21) 


(2.3.22) 


(2.3.23) 


This  set  provides  the  minimum  requirements  of  any  second-order  closure.  The  combination  of  terms  in  each 
of  these  equations  provides  a  destruction  term  for  each  correlation  which  allows  an  equilibrium  value  to 
be  reached  in  a  time  long  compared  with  A/q,  and  a  diffusion  term  which  prevents  any  excessively  sharp 
gradients  from  occurring.  Rather  than  use  models  with  a  large  nunber  of  coefficients  that  are  finely 
adjusted  to  fit  a  few  particular  flows,  we  choose  to  work  with  a  relatively  small  number.  The  appropriate 
values  of  the  coefficients  as  evaluated  by  Lewellen  (1977)  are  b  *  0.125,  A  *  0.75,  s  *  1.8,  and  vc  -  0.3. 
Results  from  other  variations  of  the  second-order  closure  set  will  be  referred  to  in  relation  to 
particular  applications. 

Scale  Determination.  To  complete  closure,  it  is  necessary  to  provide  some  means  for  determining  the 
turbulent  scales  a,  Tq,  and  involved  in  the  modeled  terms.  This  is  approached  in  different  ways  by 
various  investigators.  It  may  be  specified  empirically  based  on  the  gross  features  of  the  particular  flow 
geometry,  or  it  may  be  predicted  from  a  semlempirically  modeled,  dynamic  differential  equation.  Each  of 
these  approaches  has  some  advantages  and  disadvantages. 

The  macroscale  A  is  defined  by  Eqs.  (2.3.15  -  2.3.17).  It  is  expected  to  be  related  to  the  integral 
scale  but,  owing  to  our  choice  of  normalization,  not  equal  to  it.  It  is  also  related  to  the  mixing  length 
used  in  first  order  closure  models.  As  such,  there  is  empirical  information  that  can  be  used  in  our 
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determination  of  this  model  parameter.  It  appears  fairly  clear  that  the  scale  cannot  exceed  some  fraction 
of  the  total  spread  of  the  region  of  turbulence  and  that.  In  some  neighborhood  of  the  mall,  it  should  be 
proportional  to  the  distance  from  the  wall.  These  two  simple  Ideas,  together  with  empirical  Information 
used  to  determine  the  two  implied  constants  of  proportionality,  are  sufficient  to  permit  the  system  to 
close  with  relatively  good  nunerlcal  results  for  many  problems  (e.g.,  Donaldson,  1973;  Lewellen,  et  al., 
1974).  Others  (e.g..  Shir,  1973;  Myngaard  and  Cot6,  1974;  Mellor  and  Yamada,  1974;  Burk,  1977;  Andre, 
et  al.,  1978;  and  Brost  and  Uyngaard,  1978)  have  specified  a  completely  empirically  determined 
distribution  of  A  across  the  region  of  Interest  in  the  same  manner  as  is  done  for  first-order  mixing 
length  approaches. 

In  an  attempt  to  remove  some  of  the  arbitrariness  of  the  specification  of  A  for  different  flows,  a 
mwber  of  Investigators  have  developed  a  modeled  dynamic  equation  for  A  or  Its  equivalent.  The  starting 
points  for  such  attempts  have  varied  widely.  A  two-polnt-veloclty-correlatlon  equation  forms  the  basis 
for  the  work  of  Rotta  (1972),  Naot,  et  al.  (19731.  Donaldson  (1972),  and  Rodl  (1972).  By  forming  an 
equation  for  the  two-point  velocity  correlation  u^(x)u{(x+r)  and  integrating  to  form  an  integral  scale.  It 
appears  appropriate  to  take 


d2A 


1 


ui(x)  ui(x+  r) 


dV 

r2- 


(2.3.24) 


The  difficulty  with  this  approach  Is  that  none  of  the  terms  can  be  integrated  exactly.  All  of  the  terms 
must  be  modeled.  This  Is  also  true  if  one  starts  with  the  equation  for  the  dissipation  c,  as  Is  favored 
by  many  (e.g.,  Lumley  and  Khajeh-Nouri ,  1974;  Harlow  and  Nakayama,  1967;  and  hanjallc  and  Launder,  1972) 
or  with  the  equation  of  vorticity  fluctuations  as  suggested  by  Daly  and  Harlow  (1970)  or  Wilcox  and  Alber 
(1972).  With  the  model  chosen  in  the  section  on  Modeled  Equations,  e  »  bq3/A,  Also  the  vorticity 
fluctuations  can  be  taken  as  proportional  to  q/A.  Thus,  with  the  aid  of  the  energy  equation,  any  of  these 
approaches  may  be  reduced  to  an  equation  for  A.  As  Bradshaw  (1972),  and  Mellor  and  Herring  (1973)  have 
pointed  out,  all  of  the  resulting  A  equations  have  the  same  form; 


DA 

Ot 


a _ au. 

Si  —  ujuj  —  -  S2bq  ♦  buoyancy  terms 
q  ax-. 


+  diffusion  terms 


(2.3.25) 


The  major  difference  In  the  various  expressions  lies  In  the  construction  of  the  turbulent  diffusion  terms. 
Unfortunately,  these  turn  out  to  be  more  Important  In  the  scale  equation  than  In  the  Reynolds  stress 
equation.  The  scale  equation  as  proposed  by  Lewellen  (1977)  for  use  with  the  set  of  equations  given  In 
2.3.18  to  2.3.23  Is 
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dgiu-je  (2.3.26) 


It  Is  immediately  obvious  that  the  scale  equation  contains  much  more  arbitrariness  than  the  Reynolds 
stress  equations,  where  many  of  the  terms  were  determined  precisely  without  recourse  to  modeling  or 
coefficients.  With  such  a  large  nunber  of  coefficients  in  the  scale  equation,  a  correspondingly  large 
nuntoer  of  different  experiments  must  be  matched  concurrently  if  the  resulting  coefficients  are  to  have  any 
Invariant  validity. 

There  Is  also  the  question  as  to  what  extent  it  is  really  appropriate  to  have  the  Integral  quantity  A 
determined  by  point  values  of  the  other  variables.  Spatial  boundary  conditions  may  be  expected  to  play  a 
much  more  Important  role  in  the  determination  of  A  than  they  do  for  the  Reynolds  stress. 

Finally,  In  some  applications  It  will  be  appropriate  to  take  A9  *  A<.  *  A  while  In  other  applications 
this  will  not  be  true.  A  formulation  for  the  scalar-covariance  decay  rate,  equivalent  to  a  A,.e  equation, 
is  given  by  Lumley  (1979). 

The  model  given  in  the  preceding  sections  Is  quite  similar  to  Mellor  and  Yamada's  level  4  model 
(1974).  In  most  of  their  subsequent  applications  they  have  favored  a  condensed  version  (level  2.5)  which 
assumes  that 


(uw,  vw)  ■  -  9*  Sm  "J".  ‘J"  (2.3.27) 


we  *  -  qt  Sh  ~  (2.3.28) 


The  contracted  form  of  Eq.  (2.3.18)  Is  carried  to  determine  q2,  and  the  functions  and  Su  are  obtained 
by  assuming  a  high  Reynolds  number  local  equilibrium  In  the  remaining  components  of  the  Reynolds  stress 
and  heat  flux  equations,  i.e.,  all  gradients  of  the  turbulent  correlations  are  dropped.  This 
approximation  should  be  valid  whenever  the  time  scale  of  the  turbulence,  A/q,  1$  less  than  the 
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characteristic  tine  scale  of  the  mean  flow.  The  resulting  algebraic  relationships  between  the 
second-order  correlations  and  the  gradients  of  the  mean  variables  can  then  be  manipulated  to  give 


and 


3A,  (ww-c,q2)/q2  ♦  9A,  A,  igg  wa/q3 
m  1  +  9Aj  A2i2gp(ae/3z)/q2 


(2.3.29) 
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(Alt  BltA2.  B2.C1)  *  (0.92,  16.6.0.74,  10.1,0.08) 


Note  that  Eqs.  (2.3.29  and  2.3.30)  are  those  given  by  Mellor  and  Yamada  (1978).  When  the  same  procedure 
Is  carried  out  on  Eqs.  (2.3.18  -  2.3.23)  the  equivalent  forms  are 
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Plots  for  some  other  choices  of  Sm  and  Sjj  are  given  by  Blackadar  (1979). 


(2.3.31) 


(2.3.32) 
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III.  Similarity  Relationships 


3.1.  Mon in-Obukhov  Scaling  in  the  Surface  Layer 

The  best  understood  region  of  the  atmospheric  boundary  layer  is  the  layer  next  to  the  surface  where 
the  turbulent  fluxes  of  mass,  momentun,  and  energy  are  essentially  constant.  Rigorously  this  layer  is 
bounded  by  the  tops  of  any  surface  roughness  features  at  low  z  and  by  an  upper  height  that  remains  much 
less  than  either  the  ratio  of  the  surface  flux  of  momentum  to  the  sun  of  all  but  the  last  term  of 
Eq.  (2.1.7),  or  the  ratio  of  the  surface  heat  flux  to  the  sun  of  all  but  the  last  term  of  Eq»_(y!.1.8Ju. 
Within  this  layer  the  mean  momentum  and  energy  equations  may  be  replaced  by  the  statement  that  ujw  and  ut> 
remain  constant  equal  to  their  respective  surface  values  which  may  be  represented  as  -ui  and  u*T*.  When 
this  characteristic  velocity,  u*,  and  characteristic  temperature,  T*.  are  used  to  normalize  all  velocities 
and  temperatures  in  the  steady,  homogeneous  version  of  equations  (2.3.18-2.3.20)  the  equations  are  found 
to  be  a  function  of  only  one  variable  provided  the  turbulent  length  scale  is  a  function  of  the  same 
variable  and  the  modeled  terms  Introduce  no  new  parameters.  This  single  variable  may  be  written  as 
n  *  zxgWTpUif  5  z/L.  Physically,  the  variable  (-z/L)  is  equal  to  the  ratio  of  the  buoyant  production 
of  turbulent  kinetic  energy  to  the  production  by  shear.  The  characteristic  length  l  is  known  as  the 
Monin-Obukhov  length  following  the  work  of  Monin  and  Obukhov  (1953).  The  resulting  set  of  ordinary 
differential  equations  as  modeled  by  Lewellen  and  Teske  (1973)  may  be  written  as: 


(3.1.1) 

(3.1.2) 

(3.1.3) 

(3.1.4) 

(3.1.5) 

(3.1.6) 

(3.1.7) 
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Figure  3.1.1  -  Normalized  atmospheric  surface  layer 
gradients  as  a  function  of  stability.  Model 
calculation  from  Lewellen  and  Teske  (19/3).  Data 
from  Businger,  et  al.  (1971). 
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Figure  3.1.2  -  Vertical  velocity 
fluctuation  in  the  atmospheric 
surface  layer.  Model  calculation 
from  Lewellen  and  Teske  (1973). 
Data  from  Wyngaard,  et  al.  (1971). 
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Figure  3.1.3  -  Temperature  fluctuation  in  the 
atmospheric  surface  layer.  Model  calculation  from 
Lewellen  and  Teske  (1973).  Data  from  Wyngaard,  et 
al.  (1971). 


These  equations  may  be  integrated  to  yield  the  _igven  normalized  variables:  »  =  (rz/u*)(aU/az), 

vh  =  (K2/T*)(ae/3z),  ow/u*.  ou/u*,  ov/u*.  uu/u*T*.  and  o„/T*.  The  resulting  functions  for 
some  of  these  variables  are  compared  with  field  observations  in  Figures  3.1.1  -  3.1.3.  In  these 
computations  the  turbulent  scale  A  was  set  equal  to  0.65z,  except  it  was  not  allowed  to  exceed  an  upper 
bound  of  0.2L  in  stable  flow.  This  particular  bound  on  A  was  chosen  to  yield  the  observed  critical 
Richardson  nunber  of  0.21.  This  bound  was  also  used  in  determining  the  coefficient  of  the  buoyant  term  in 
the  scale  equation,  Eq.  (2.3.26). 

The  normalized  velocity  and  temperature  gradients.  Figure  3.1.1,  show  very  good  agreement  between  the 
model  calculations  and  the  data.  A  positive  surface  heat  flux  (z/L  <  0)  leads  to  increased  turbulence  and 
decreased  gradients,  while  a  negative  surface  heat  flux  (z/L  >  0)  leads  to  a  damping  of  the  turbulence  and 
increased  gradients.  Under  neutral  conditions  (z/L  *  0),  the  ratio  of  the  effective  diffusivity  for 
momentum  to  that  for  temperature  is  3/4.  In  the  model,  this  choice  was  determined  by  the  value  chosen  for 
the  model  constant  A  in  Eq.  (2.3.12). 

The  vertical  velocity  and  temperature  fluctuations  are  compared  with  the  data  of  Wyngaard, 
et  al.  (1971),  in  Figures  3.1.2  and  3.1.3.  The  agreement  for  the  vertical  velocity  is  very  good  but  that 
for  the  temperature  fluctuations  leaves  a  bit  to  bg  desired.  As  pointed  out  by  Wyngaard,  et  al.  (1971), 
there  is  considerable  uncertainty  in  (U7) l/zu*/wu  at  z/L  «  U  because  both  variables  go  to  zero. 

Mellor  (1973)  made  a  somewhat  similar  calculation  of  these  surface-layer  similarity  functions.  The 
major  difference  is  that  Mellor  eliminated  the  diffusion  terms  so  that  the  differential  equations  reduce 
to  an  algebraic  set.  In  this  case  the  length  scale  may  be  normalized  out  of  the  problem,  eliminating  any 
need  for,  or  any  possibility  of,  incorporating  an  influence  of  stratification  on  the  scale.  His  model 
coefficients  must  then  be  chosen  so  that  they  will  match  the  critical  Richardson  nunber. 

Although  this  surface  layer  similarity  also  leads  to  variations  of  ou  and  ov  as  a  function  of  z/L, 
observations  of  these  variables  do  not  correlate  with  z/L.  These  horizontal  wind  variances  are  more 
dependent  on  variables  which  control  the  large  eddy  structure  at  larger  heights  above  the  boundary  layer. 
This  behavior  can  be  incorporated  within  the  surface  layer  similarity  if  the  modeled  higher  order 
correlations  in  Equations  (2.3.18-2.3.20)  are  allowed  to  be  a  function  of  two  scales,  rather  than  one.  A 
model  which  allows  the  horizontal  wind  variance  to  be  more  a  function  of  boundary  layer  thickness  than  it 
is  of  height  z  has  been  developed  by  Lewellen  and  Sandri  (1980).  The  results  of  this  model  are  compared 
with  observations  in  the  unstable  layer  (see  Section  3.3)  where  the  horizontal  variance  has  been  found  to 
correlate  with  z-j /L . 

Although  the  normalized  velocity  and  temperature  gradients  in  the  surface  layer  can  be  completely 
determined  by  z/L  only,  the  velocity  and  temperature  at  any  height  will  also  be  a  function  of  surface 
roughness.  Useful  analytic  functions  may  be  obtained  by  analytically  fitting  the  empirical  curves  given 
in  Figures  3.1.1  and  3.1.2  [e.g.,  Buslnger,  et  al.  (1971) j  and  integrating  as  carried  out  by 
Paulson  (1970)  the  resulting  expressions  for  U  and  y  are: 
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The  species  concentration  may  be  expected  to  have  a  variation  similar  to  the  temperature. 

Many  variations  on  the  formulation  of  scaling  In  the  surface  layer  have  been  given  in  the  literature, 
e.g.,  Blackadar,  1979;  Bienkowski,  1979;  Herbert  and  Panhaus,  1979;  and,  Panhaus  and  Herbert.  1979. 


3.2.  Neutral,  Steady-State,  PBt  Scaling 

The  governing  equations  as  given  in  Chapter  II  reduce  in  this  case  to  the  two  components  of  the 
momentum  equation: 


fv  -  (uw) 

(3.2.1) 

f(Ug-U)  -  ^  (VW) 

(3.2.2) 

and  the  six  components  of  the  Reynolds  stress  equation.  Here  the  coordinate  system  has  been  aligned  with 
the  geostrophic  wind  so  that  ap/»x  *  0  and  *P/3y  *  -fUg.  The  appropriate  boundary  conditions  for  this 
case  are  to  have  the  mean  wind  approach  the  geostrophic  wind,  Ug,  as  z  *  -  and  to  go  to  zero  at  the 
effective  aerodynamic  roughness  height,  z0,  while  all  the  turbulent  quantities  go  to  zero  as  z  *  -  and 
approach  the  surface  layer  similarity  values  given  in  Section  3.1  as  z  ♦  z0. 


Under  these  conditions,  the  governing  equations  can  be  normalized  so  that  only  one  parameter,  the 
Rossby  nunber  Ro  «  Ug/z0f,  is  needed  to  determine  the  flow.  The  resulting  relative  variation  with  respect 
to  Rossby  number  shows  how  the  flow  responds  to  changes  in  geostrophic  wind,  surface  roughness,  and 
latitude  which  determines  f.  Figures  3.2.1  and  3.2.2  show  the  mean  wind  profiles  normalized  by  the 
geostrophic  wind  as  a  function  of  the  normalized  height  for  different  values  of  Ro  as  given  by  a  numerical 
solution  to  the  Reynolds  stress  equations  as  modeled  in  Eq.  (2.3.18)  (Lewellen  and  Williamson,  1976).  The 
change  in  wind  direction  with  height  is  less  than  it  would  be  in  an  Ekman  Spiral,  which  is  the  solution  of 
the  analogous  laminar  flow  problem.  The  vertical  velocity  variance  is  shown  in  Figure  3.2.3  for  the  same 
conditions. 

The  vertical  coordinate  for  these  plots  Is  normalized  in  such  a  way  that  for  typical  values  of 
f  »  10-Sec-1  and  Ug  *  10  m/sec,  z  may  be  read  directly  in  Kms.  For  other  values  of  Ug/f,  the  plot  must 
be  read  as  a  dimensionless  plot  and  the  appropriate  physical  z  determined.  An  Increase  Tn  the  geostrophic 
wind  will  Increase  Ro  and  reduce  the  dimensionless  boundary  layer  thickness  but  the  physical  boundary 
layer  thickness  will  Increase  since  the  effect  in  the  coordinate  transformation  is  much  stronger  than  that 
produced  by  an  equal  change  in  Ro.  Likewise  the  boundary  layer  will  become  thicker  as  the  equator  is 
approached,  for  equal  values  of  Ug/z0.  Increases  in  the  surface  roughness  always  increase  the  boundary 
layer  thickness  in  this  limit. 


Rossby  Number  Variation  -  U 


Rossby  Number  Variation  -  V 


Figure  3.2.1  -  Influence  of  Rossby 
number  on  the  mean  wind  component  in 
the  direction  of  the  geostrophic 
wind  (normalized  altitude  converts 
to  z  in  Km  for  f  *  lO'Sec'1  and 
Ug  «  10  m/sec) . 


Figure  3.2.2  -  influence  of  Rossby 
number  on  the  mean  wind  component  in 
the  direction  normal  to  the 
geostrophic  wind  for  the  conditions 
of  Figure  3.2.1. 
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Figure  3.2.4  -  Surface  shear  stress 
as  a  function  of  Rossby  number  for 
the  neutral  planetary  boundary 
layer. 


Figure  LI.*?.  "  Inf1uence  of  Rossl*  on  the  variance  of  the  vertical  velocity  for  the  conditions  of 
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Figure  3.2.5.  -  Angle  between  the 
surface  wind  and  the  geostrophlc 
wind  as  a  function  of  Rossby 
number  for  the  neutral  planetary 
boundary  layer 


Figure  3.2.6  -  Comparison  of  the 
solution  predicted  by  Eqs.  (3.2.3 
and  3.2.4)  (— )  with  that  given  for 
the  complete  second-order  closure 
model  in  Figures  3.2.1  and  3.2.2. 
( - )  U0  *  0-1.  Ug/fZo  *  l°b)- 


The  surface  drag  coefficient,  defined  as  Ul/U^,  is  a  function  only  of  Ro  in  this  limiting  case.  It 
is  given  in  Figure  3.2.4  as  predicted  by  the  present  model.  A  summary  of  empirical  correlations  for  this 
variable  is  given  by  Me Bean,  et  a 1.  (1979).  The  data  show  considerable  scatter  since  the  conditions 
necessary  for  this  steady  state,  neutral  profile  are  seldom  completely  satisfied  in  nature.  The  model 
result  shown  falls  within  the  range  of  the  empirical  correlations.  The  slight  decrease  in  drag 
coefficient  with  increasing  Ro  is  consistent  with  a  decreasing  ratio  of  characteristic  turbulent 
macroscale  A  to  boundary  layer  thickness.  As  shown  in  Figure  3.2.5,  the  angle  between  the  surface  wind 
and  the  geostrophlc  wind  also  decreases  with  increasing  Ro. 

The  logarithmic  variation  for  wind  speeds  in  the  surface  layer  given  by  Eq.  (3.1. a)  under  neutral 
conditions  is  valid  to  a  much  higher  altitude  for  U  than  it  is  for  V.  The  departure  from  tn  z/zg  is 
forced  by  the  geostrophic  pressure  gradient  fUg  acting  parallel  to  v.  The  logarithmic  term  is  the  first 
term  in  a  consistent  expansion  for  small  z.  When  the  expansion  is  continued  to  bring  in  the  influence  of 
the  pressure  gradient  the  two  velocity  components  may  be  written  as  (lewellen,  1977): 
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(3.2.3) 
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Equations  (3.2.3)  and  (3.2.4)  agree  with  Figures  3.2.1  and  3.2.2  to  heights  well  above  that  at  which 
occurs,  as  seen  in  Figure  3.2.6. 


3.3.  Convectively  Mixed  Layer  Scaling 

The  unstable  boundary  layer  which  is  dominated  by  buoyantly  driven  turbulence  has  proven  to  be 
amenable  to  a  similarity  description.  In  a  series  of  publications  involving  analysis,  nunerical  studies 
and  experimental  observations  Deardorff  and  Hi  Ills  (Deardorff,  1970a, b,  1972,  1974a, b;  and,  Willis  and 
Deardorff,  1974)  have  shown  that  the  mean  temperature,  and  several  of  the  turbulent  correlations  are 
functions  only  of  z/z\  and  of  the  characteristic  convective  velocity  W*.  This  characteristic  velocity  is 
defined  as 
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U* 


(3.3.1) 


where  w«0  is  the  heat  flux  near  the  lower  surface  and  Zj  is  the  depth  of  the  mixed  layer.  Such  a 
similarity  holds  in  the  unstable  PBL  above  the  height  -L,  where  L  is  the  Monin-Obukhov  length  discussed  in 
Section  3.1. 

The  similarity  functions  for  this  buoyantly  driven  turbulence  may  be  obtained  from 
Eqs.  (2.3.18  -  2.3.20)  by  setting  the  mean  velocity  equal  to  zero,  imposing  a  fixed  surface  heat  flux  at 
the  bottom,  and  zero  turbulence  and  a  stable  lapse  rate  at  the  top  of  the  domain.  The  equations  for  this 
limit  are 


(3.3.2) 

(3.3.3) 

(3.3.4) 

(3.3.5) 

(3.3.6) 

(3.3.7) 


The  resulting  similarity  functions  as  given  by  Lewellen,  Teske  and  Donaldson  (1976)  are  plotted  in 
Figures  3.3.1  -  3.3.3. 

Figure  3,3.1  compares  the  model  predictions  for  the  normalized  vertical  velocity  fluctuations  with 
the  experimental  observations.  The  model  predictions  assume  no  mean  velocity  and  a  constant  surface  heat 
flux  to  simulate  the  laboratory  experiment.  The  agreement  between  predictions  and  observations  is 
heartening,  particularly  since  no  model  constant  adjustments  were  involved. 


Figure  3.3.1  -  Similarity  profiles  of  vertical 
velocity  fluctuation  In  an  unstable  mixed  layer; 
model  calculation  from  Lewellen,  Teske,  and 
Donaldson  (1976);  data  from  Willis  and 
Oeardorff  (1974). 


Figure  3.3.2  -  Similarity  profiles 
of  horizontal  velocity  fluctuations 
for  the  conditions  of  Figure  3.3.1. 
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The  horizontal  velocity  fluctuations  are  given  in  Figure  3.3.2.  Here  the  comparison  with  experiment 
not  as  close  since  the  observations  show  a  peak  near  the  lower  surface  that  is  not  evidenced  in  the 
model  results.  The  distributions  of  temperature  fluctuations  are  shown  in  Figure  3.3.3.  The  agreement 
between  model  and  experiment  is  very  good  except  near  the  top  of  the  mixed  layer  where  the  observations 
show  a  much  stronger  local  maximum  than  is  predicted.  This  is,  at  least  in  part,  due  to  the  existence  of 
internal  gravity  waves  in  this  region. 

Favorable  comparisons  between  model  results  and  observations  exist  in  Figures  3.3.1  -  3.3.3  in  spite 
of  the  fact  that  the  modeled  vertical  flux  of  kinetic  energy  is  significantly  different  from  the 
experimental  observations,  as  shown  in  Figure  3.3.4.  This  is  consistent  with  our  basic  assumption  that 
the  second-order  correlations  calculated  from  their  conservation  equations  are  given  more  accurately  than 
the  accuracy  to  which  the  third-order  correlations  are  approximated. 


Figure  3.3.3  -  Similarity  profiles  Figure  3.3.4  -  Similarity  profiles 
of  temperature  fluctuations  for  the  of  vertical  flux  of  kinetic  energy 
conditions  of  Figure  3.3.1.  for  the  conditions  of  Figure  3.3.1. 


As  shown  in  Figure  3.3.5,  Zeman  and  Lumley  (1976)  obtain  much  better  agreement  between  predictions  of 
the  flux  of  kinetic  energy  and  the  observations  of  this  third-order  quantity,  using  their  second-order 
model  designed  specifically  for  buoyancy  driven  mixed  layers.  However,  the  comparisons  for  the 
second-order  correlations,  shown  in  Figures  3.3.6  -  3.3.7  are  not  substantially  different  from  that 
obtained  with  the  simpler  model.  The  significant  discrepancy  between  model  results  and  observations  in 
the  horizontal  fluctuations  near  the  surface  is  unaffected  by  the  more  rigorous  modeling  of  the 
third-order  terms.  Similar  results  for  the  second-order  correlations  of  turbulence  in  the  convectively 
mixed  layer  have  also  been  obtained  by  Mellor  and  Vamada  (1976),  and  Sun  and  Ogura  (1960). 


Figure  3.3.5  -  Zeman  and  Lumley's  ,  ,  ,, 

(1976)  predicted  fluxes  of  turbulent  Figure  3.3.6  -  Vertical  velocity  fluctuations 
energy  compared  with  data:  solid  as  calculated  by  Zeman  and  Lumley  (1976) 
circles  and  triangles,  Willis  and  compared  with  data:  solid  circles  and 
Deardorff  (1974);  open  circles  and  triangles,  Willis  and  Deardorff  (1974);  open 
triangles,  Lenschow  (1970,  1974)  and  circles  and  triangles,  Lenschow  (1970,  1974)  and 
Lenschow  and  Johnson  (1968);  and  Lenschow  and  Johnson  (I960);  and  squares, 
squares,  Telford  and  Warner  (1964),  Tel  ford  and  Warner  (1964) . 
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Figure  3.3.7  -  Zeman  and  Lumley's 
(1976)  predicted  temperature 
variance  compared  with  data:  solid 
circles.  Case  SI  Willis  and 
Deardorff  (1974);  open  circles  arid 
triangles,  Lenschow  (1970,  1974); 
squares,  Telford  and  Warner  (1964). 


Figure  3.3.8  -  Evolution  of  the 
velocity  spectra  with  height  as 
given  by  a  composite  plot  of  a 
number  of  experiments  (Kaimal, 
1978). 


The  principal  difficulty  with  modeling  the  horizontal  fluctuations  appears  to  be  the  anisotropic 
nature  of  the  turbulence  scales,  as  mentioned  in  Section  3.1.  This  is  demonstrated  in  the  composite  plot 
of  the  velocity  spectra  as  a  function  of  height  given  in  Figure  3.3.8  from  Kaimal  (1978).  The  vertical 
velocity  component  clearly  shows  the  scale  of  the  eddies  with  the  highest  energy  increasing  in  length 
proportional  to  increases  in  height,  as  assumed  in  Section  3.1.  However,  the  corresponding  scale  of  the 
horizontal  velocity  fluctuations  shows  almost  no  change  with  height.  Instead,  tn*  he  izontal  velocity 
fluctuations  appear  to  scale  with  the  mixed-layer  height,  as  exhibited  by  Wyngaaru  (1979'. 


Figure  3.3.9  -  Dimensionless  rms  horizontal  velocity  fluctuation  as  a  function  of  the 
ratio  of  the  mixed  layer  height  to  the  Monin-Qbukhov  length.  Two  scale  model  calculation 
from  Lewel len  and  Sandri  (1980);  data  from  Wyngaard  (1979). 


Lewellen  and  Sandri  (1980)  have  attempted  a  modification  to  the  model  presented  in  Section  2.3  which 
permits  the  horizontal  macro-length  scale  to  be  different  from  the  vertical  macro-length  scale.  The 
vertical  scale  is  identified  with  the  previous  single  scale  as  modeled  in  Section  2.3  and  the  horizontal 
scale  is  allowed  to  vary  with  the  mixed- layer  height. 

The  results  of  the  2-scale  model  for  the  horizontal  velocity  variance  are  contrasted  with  those  of 
the  single-scale  model  in  Figure  3.3.9.  The  single-scale  model  predicts  a  variation  with  respect  to  z/L 
only,  which  is  not  observed  in  the  data.  Rather,  the  data  are  found  to  be  relatively  independent  of  z/L 
and  to  depend  most  strongly  on  z^/L.  The  results  of  the  2-scale  model  are  consistent  with  the  data.  The 
challenge  of  providing  a  consistent  forwjlation  for  a  general  2-scale  model  remains  to  be  solved. 

3.4.  The  Stable  Ekman  Boundary  Layer 

When  the  surface  is  cooler  than  the  air  above  it  a  positive  (stable)  temperature  gradient  is 
typically  established.  This  has  a  strong  dynamical  effect  on  the  turbulent  structure  of  the  layer.  The 
most  important  parameter  governing  these  effects  is  a  Richardson  ntxnber  which  measures  the  ratio  of 
buoyancy  and  inertia  forces.  Fundamental  studies  of  shear  flow  in  a  stably  stratified  mediun  (e.g.. 
Miles,  1961;  and  Howard,  1961)  show  that  when 
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linear  disturbances  are  damped  out.  When  Ri  is  less  than  this  critical  value,  turbulence  may  be  expected 
to  grow  driven  by  the  shear  and  Ri  be  increased  by  the  decreases  in  gradients  forced  by  the  turbulence. 
When  Ri  >  Ricrit*  turbulence  may  be  expected  to  be  damped  by  buoyancy  which  in  the  presence  of  any  heat  or 
momentum  transport  will  lead  to  increased  gradient  and  decreases  in  Ri.  Thus,  as  suggested  by  Long  (1976) 
natural  stratified  shear  layers  should  tend  towards  an  equilibrium  at  Ricrit*  average  values  of  the 
gradients  are  used,  then 
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and  setting  this  equal  to  Ricrit  leads  to 
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consistent  with  the  estimate  of  layer  thickness  given  in  Eq.  (2.2.9).  In  general,  Eq.  (3.4.3)  provides 
for  an  evolving  boundary  layer  with  and  6  varying. 

When  Ri4  ■  R’crit  and  6  satisfies  the  momentum  equation  so  that  6  ~  q/f,  then  a  steady  equilibrium 
value  of  a  is  possible.  The  thickness  of  this  equilibrium  Ekman  layer  may  be  given  following 
Zil itinkevich  (1972)  as 


6  «  const  (u*L/f) 


(3.4.4) 


It  may  be  shown  to  be  consistent  with  that  given  by  Eq.  (2.2.8). 

The  approach  to  a  stable  equilibrium,  as  given  by  Wyngaard  (1975),  using  his  version  of  second-order 
closure,  is  given  in  Figure  3.4.1  as  a  function  of  cooling  rate.  Three  to  six  hours  is  the  time  required 
to  reach  the  equilibrium  state  according  to  these  calculations.  This  long  transition  time  implies  that 
such  equilibrium  stable  layers  will  seldom  be  achieved  because  boundary  conditions  seldom  remain  constant 
this  long  in  the  atmosphere.  Nonetheless,  it  is  instructive  to  look  at  the  wind  and  turbulence  profiles 
expected  under  such  steady  conditions.  Under  such  conditions  the  normalized  profiles  of  wind, 
temperature,  and  turbulent  correlations  should  be  a  function  of  only  two  parameters,  the  Rossby  number, 
Ug/Zof.  and  a  stability  parameter.  The  most  frequent  choice  of  the  stability  parameter  is  u  *  u*/fL 
(ZIHtinkevich,  1972).  In  terms  of  this  parameter,  Eq.  (3.4.4)  may  be  written  as 
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Near  similarity  can  be  achieved  using  Lul/i  as  a  characteristic  length  to  normalize  the  governing 
equations  for  the  equilibrium  case  of  strong  stability.  When  u*  is  used  to  normalize  all  turbulent 
velocities,  e»  to  normalize  all  temperatures,  and  u»u*'^  to  normalize  the  mean  velocity  the  normalized 
Reynold's  stress  equation  may  be  written  as 
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It  is  possible  to  show  from  the  scale  equation  that  a  should  be  proportional  to  L  under  strongly  stable 
conditions  so  that  X  ~  tT1/2.  Thus  the  diffusion  term  in  Eq.  (3.4.6)  is  the  only  term  that  keeps  the 
Reynolds  stress  equation  from  becoming  independent  of  u,  and  Ro.  For  large  u  the  diffusion  term  should  be 
negligible  since  it  is  order  u"1  in  comparison  with  the  other  terms.  Likewise  the  momentum  equations  can 
be  made  independent  of  u  and  Ro.  The  derivative  of  the  normalized  equations  may  be  written  as 
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The  appropriate  boundary  conditions  for  Eqs.  (3.4.6  -  3.4.8)  are  also  independent  of  m  and  Ro. 
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(3.4.9) 


The  only  remaining  requirement  for  similarity  is  that  the  thermal  boundary  conditions  be  specifiable  in  a 
manner  which  leaves  the  heat  flux  equation  similar. 


The  only  difficulty  lies  with  the  energy  equation 
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(3.4.10) 
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which  cannot  be  independent  of  time  and  still  satisfy  boundary  conditions  which  permit  the  heat  flux  to 
have  a  zero  value  at  1  *  ».  However,  this  can  nearly  be  achieved  by  postulating  a  fixed  cooling  rate 
which  allows  the  heat  flux  to  vary  linearly  with  ?  across  the  bulk  of  the  boundary  layer. 

The  resulting  distributions  for  the  Reynolds  stresses  for  this  “equilibrium"  stable  layer  are  given 
In  Figure  3.4.2  as  calculated  by  Wyngaard  (1975).  The  thickness  of  the  boundary  layer,  i.e.,  the  constant 
In  Eq.  (3.4.5),  Is  fairly  sensitive  to  the  details  of  the  model.  Wyngaard  obtained  a  value  of 
approximately  0.2.  The  corresponding  values  of  this  constant  from  other  models  are  0.55  (Businger  and 
Arya,  1974),  0.4  (Brost  and  Wyngaard,  1978),  0.27  (Rao  and  Snodgrass,  1979),  and  approximately  (0.5)  from 
the  model  given  in  Eqs.  (3.4.o  -  3.4.12).  Attempts  to  determine  this  mmber  from  observations  have  not 
been  successful,  presumably  due  to  the  scatter  caused  by  such  phenomena  as  advective  effects,  radiative 
flux  divergence,  barocl inicity,  or  unsteady  influences  to  be  discussed  in  Chapter  IV.  There  is  also  some 
arbitrariness  in  the  definition  of  6.  For  the  preceeding  analysis  4  has  been  defined  as  the  height  at 
which  the  shear  stress  falls  to  5%  of  Its  surface  value.  This  height  is  well  below  the  height  at  which 
the  velocity  goes  to  the  geostrophlc  value  due  to  its  oscillatory  approach  to  «. 
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Figure  3.4.1  -  The  behavior  of  u«  in  an  evolving 
nocturnal  P8L,  as  predicted  by  the  numerical  model  of 
Wyngaard  (1975).  Various  cooling  rates  of  the 
near- surface  air  were  used,  as  shown. 
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Figure  3.4.2  -  The  steady-state  turbulence  profiles 
in  the  nocturnal  PBL,  as  found  from  the  model 
calculations  of  Wyngaard  (1975)  (Northern 
Hemisphere). 


Although  the  turbulence  correlations  can  be  presented  in  a  normalized  plot  independent  of  u  and  Ro, 
the  wind  profiles  themselves  cannot.  They  remain  a  function  of  both  Rossby  number  and  stability.  Wind 
profiles  for  particular  values  of  these  parameters  are  given  in  Section  3.5. 

It  must  be  recognized  that  this  treatment  of  the  stable  layer  is  based  on  a  gross  simplification  of 
its  true  dynamics,  since  it  ignores  any  role  of  Internal  gravity  waves.  Recent  reviews  of  Internal 
gravity  waves  in  the  stable  boundary  layer  have  been  given  by  Chimonas  (1980)  and  Einaudi  (1980).  Such 
waves  can  transport  momentum  vertically  across  the  layer  or  horizontally  within  the  layer  and  Interact 
strongly  with  the  turbulence.  The  model  presented  In  Chapter  11  is  based  on  ensemble  averaging  over 
fluctuations  In  the  flow.  Equations  (2.1.10  -  2.1.12)  are  as  valid  for  wave  fluctuations  as  they  are  for 
turbulent  fluctuations.  However,  the  modeled  terms  in  Eqs.  (2.3.18  -  2.3.23)  have  been  specifically 
wdeled  for  fully  turbulent  flow  and  thus  cannot  be  expected  to  be  valid  for  coherent  wave  fluctuations. 
This  remains  a  challenge  for  the  future  to  generalize  the  modeling  to  represent  both  wave  and  turbulent 
fluctuations. 
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3.5.  Stability  Variations 

Although  we  expect  the  dynamics  of  the  changing  boundary  conditions  to  be  important  for  most  real 
stability  variations  in  the  atmosphere,  we  will  first  consider  the  quasi -steady  case  where  the  surface 
heat  flux  is  held  constant.  Figures  3.5.1  and  3.5.2  show  the  wind  profiles  for  four  different  stability 
conditions  with  Ro  held  constant  at  10b. 


Figure  3.5.1  -  Influence  of 
Richardson  number  on  the 
quasi-steady  distributions  of  mean 
wind  component  in  the  direction  of 
the  geostrophic  wind  (Ro  *  10b); 
Lewellen  and  Williamson  (1976). 
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Figure  3.5.2  -  Influence  of 
Richardson  number  on  the 
quasi-steady  distributions  of  mean 
wind  component  in  the  direction 
normal  to  the  geostrophic  wind  for 
the  conditions  of  Figure  3.5.1. 


A  measure  of  the  stability  which  is  sometimes  more  convenient  than  v  is  that  of  a  bulk  Richardson 
number  based  on  the  temperature  and  velocity  differences  between  the  surface  and  an  arbitrary  height  such 
as  10  meters. 
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From  Eqs.  (3.1.10  and  3.1.11)  it  is  possible  to  see  that  the  most  stable  value  possible  for  Rilu  in  the 
surface  layer  is  (4.7 )" 1  =>  0.21,  the  critical  Richardson  number  given  by  this  model.  Larger  values  of 
Rilumean  that  the  surface  layer  is  less  than  10  meters  thick.  When  the  bulk  Richardson  number 
defined  by  Eq.  (3.5.1)  is  positive,  the  flow  is  stable  and  the  turbulence  is  damped.  As  shown  in 
Figures  3.5.1  and  3.5.2,  Ril0  *  0.0135  produces  a  boundary  layer  which  is  significantly  smaller  than  the 
neutral  case  of  Ri10  *  0,  as  should  be  expected  for  damped  turbulence.  This  also  produces  a  stronger 
transverse  wind  component,  as  seen  in  Figure  3.5.2.  Although  the  wind  shear  at  the  surface  is  reduced  by 
the  increase  in  stability,  the  wind  shear  a  short  distance  above  the  surface  is  significantly  increased. 
Thus,  from  the  point  of  view  of  the  wind  shear  which  an  aircraft  may  encounter,  an  increase  in  Rim  is 
quite  likely  to  be  associated  with  increased  wind  shear,  as  will  be  discussed  in  Section  6.2. 


As  outlined  in  Section  3.3,  the  specification  of  Rim  alone  is  not  sufficient  to  determine  the 
unstable  profile  (Rim  <  0)  even  in  the  quasi-steady  case.  Figures  3,5.1  and  3.5.2  show  two  wind  profiles 
for  nearly  the  same  negative  value  of  Rim  and  two  quite  different  values  of  the  inversion  height.  The 
unstable  layer,  driven  by  a  positive  surface  heat  flux,  must  be  capped  by  a  stable  temperature  inversion 
layer  if  the  boundary  layer  is  to  reach  a  quasi-steady  state.  The  figures  show  that  in  addition  to  the 
shear  layer  at  the  surface,  a  shear  layer  is  also  developed  In  the  vicinity  of  the  temperature  inversion. 
As  may  be  noted  in  the  discussion  of  the  boundary  conditions  in  Section  3.3,  i\  is  not  a  rigid  lid,  but  is 
defined  as  the  altitude  at  which  the  maximum  values  of  u*  occur  as  the  temperature  Inversion  is 
asymptotically  approached.  Between  the  surface  shear  layer  and  the  inversion  layer,  the  layer  is 
relatively  uniformly  mixed. 
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Figure  3.5.3  -  Vertical  velocity  variance  for 
unstable  conditions  compared  with  the  limiting 
profile  for  free  convection.  Lewellen  and  Williamson 
(1976). 

The  vertical  velocity  variance  is  shown  in  Figure  3.5.3  for  the  unstable  cases,  but  with  the 
normalization  which  is  appropriate  for  free  convection  when  there  is  no  mean  wind  shear.  The  limiting 
profile  for  R1lu  ♦  --was  shown,  in  Section  3.3,  to  agree  very  well  with  laboratory  simulations  of  this 
limit.  Figure  3.5.3  shows  that  it  does  not  require  a  very  negative  Rim  to  approach  this  limiting 
distribution.  The  corresponding  values  of  the  Zilitinkevich  parameter,  m,  are  also  given. 


3.6.  Flow  within  a  Surface  Canopy 


A  canopy  of  vegetation  presents  a  complex  lower  boundary  for  atmospheric  flows.  For  flow  well  above 
this  canopy,  it  is  usually  adequate  to  characterize  the  boundary  in  terms  of  only  an  aerodynamic 
roughness,  z0.  But  when  one  is  interested  in  the  flow  within  the  canopy  or  immediately  above  it,  a  more 
detailed  representation  is  required. 


Figure  3.6.1  -  Wind  and  turbulence  profiles  in  and  above  a  corn  canopy.  Comparison  of 
model  predictions  (Lewellen  and  Sheng,  1980)  with  data  from  Shaw,  et  al.  (1974). 
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Figure  3.6.1  -  Wind  and  turbulence  profiles  in 
and  above  a  corn  canopy.  Comparison  of  model 
predictions  (Lewellen  and  Sheng,  1980)  with  data 
from  Shaw,  et  al.  (1974). 

Second-order  closure  models  for  canopy  flow  have  recently  been  given  by  Wilson  and  Shaw  (1977),  and 
Lewellen  and  Sheng  (1980).  The  principal  difference  between  these  two  models  is  that  the  latter  consider 
heat  and  species  transport  as  well  as  momentwn  transport.  They  also  use  a  somewhat  more  general 
representation  of  the  drag  per  unit  volume  of  the  vegetation.  The  model  is  designed  to  predict  the 
variation  in  surface  layer  heat  and  species  transport  as  a  function  of  surface  Reynolds  number  (U*Zg/v), 
Prandtl  number,  K/v,  Scimidt  number,  Kc/v,  and  plant  area  density  distribution.  The  comparison  of  the 
mean  distributions  of  velocity,  Reynolds  stress  and  velocity  variance  predicted  by  Lewellen  and 
Sheng  (1980)  are  compared  with  some  limited  data  for  a  corn  canopy  (Shaw,  et  al.,  1974)  in  Figure  3.6.1. 

The  distributions  shown  in  Figure  3.6.1  are  only  appropriate  for  vegetal  canopies  with  plant  area 
distribution  similar  to  that  of  a  mature  corn  field.  If  the  plant  area  distribution  is  thicker  near  the 
top  of  the  canopy,  the  effective  displacement  height  of  the  mean  wind  will  be  increased.  While  a  less 
dense  canopy,  such  as  a  winter  deciduous  forest,  will  lead  to  larger  values  of  u*  with  a  smaller 
displacement  height.  An  approximately  similar  way  of  presenting  the  wind  profiles  for  quite  different 
canopy  structures  is  given  in  Figure  3.6.2  where  (U  -  Uy, ) /u*  is  plotted  as  a  function  of  (z  -  d)/lh  -  d). 
Both  model  calculations  and  field  data  show  approximately  similar  shapes  when  presented  this  way. 
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Figure  3.6.2  -  Similarity  profile  of  gaseous 
concentration  (1)  and  wind  velocity  (2)  predicted  by 
Lewellen  and  Sheng  (1980)  and  measured  wind  velocity 
as  given  by  Cionco  (1971). 
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More  details  on  the  flow  within  particular  canopies  are  given  in  the  two  volumes  edited  by 
Montieth  (1976)  on  vegetation  and  the  atmosphere.  Other  information  may  also  be  found  in  the  work  of 
Cionco  (1972)  and  Raupach  and  Thom  (1981). 

The  surface  layer  formulation  given  in  Section  3.1  is  generally  expected  to  hold  above  the  canopy 
with  z  measured  from  the  effective  displacement  height  within  the  canopy.  However,  as  shown  by 
Figure  3.6.3,  when  d  is  chosen  as  the  displacement  height  for  momentum,  the  inverse  gradient  functions, 
,  may  show  considerable  variation  from  the  typical  behavior.  These  data  show  effective  £  which 
ar4  double  their  expected  values  above  a  forest  under  slightly  unstable  conditions.  The  data  sn&m  were 
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Figure  3.6.3  -  Stability  dependence  of  inverse 
influence  functions  ^  r,  measured  over  Thetford 
forest  (  □)  and  Uriarr’a  ’forest  (  ■  ).  The  smooth 
curves  represent  typical  inertial-sublayer  behavior. 
Raupach  and  Thom  (1981). 


taken  approximately  5  meters  above  two  different  forests  of  height  h  x  16  m  and  a  mean  tree  spacing  of  i 
to  3  m.  This  indicates  that  the  distribution  of  sources  and  sinks  for  momentum,  heat,  and  water  vapor  may 
be  quite  different  within  a  forest  leading  to  characteristic  length  scales  for  the  fluctuations  of  heat 
and  water  vapor  different  from  those  for  the  fluctuations  of  momentum.  Thus  one  must  be  very  cautious  in 
applying  simple  scaling  laws  to  flow  within,  or  just  above,  real  canopies. 


IV.  Interaction  of  Turbulent  Transport 
with  other  Physical  Mechanisms 


4.1  Diurnal  Variations  in  the  Planetary  Boundary  Layer 

Over  homogeneous  terrain,  the  atmospheric  boundary  layer  may  be  considered  a  function  of  time  and 
altitude  only,  when  horizontal  variations  in  the  boundary  conditions  at  the  top  of  the  boundary  layer  are 
considered  negligible.  Under  such  ideal  conditions,  the  boundary  layer  may  be  expected  to  pass  through 
the  stability  conditions  considered  in  the  last  chapter  during  the  course  of  a  typical  day.  Since  thermal 
energy  is  principally  received  and  emitted  from  the  surface,  we  expect  unstable  conditions  to  dominate  in 
the  afternoon  hours  and  stable  conditions  to  dominate  during  the  late  night  hours.  The  results  of 
calculations  for  a  fixed  geostrophic  wind  and  upper  level  temperature  lapse  rate  with  a  cyclic  surface 
heat  flux-  approximating  conditions  over  the  Midwestern  plains  during  summer  were  presented  by 
Lewellen,  et  al.  (1974),  These  results  were  obtained  using  the  quasi  equilibrium  approach  outlined  in 
Section  2.3,  and  the  gross  feature  approach  to  the  scale.  Results  obtained  using  the  full  equations  were 
presented  by  Lewellen  (1977).  There  is  relatively  little  difference  in  the  nunerical  results. 


0  9  C  9  12  15  1*  91  24 

TMC  (HRS  ) 

Figure  4.1.1  -  Isopleths  of  turbulent  energy  in 
diurnal  boundary  layer.  Lewellen  (1977). 
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Figure  4.1.2  -  Vertical  distributions  of  U  taken  at 
select  times  from  a  model  simulation  of  a  midwestern 
summer  day.  Lewellen,  et  al.  (1974). 
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Figure  4.1.3  -  The  vertical  distribution  of  V  taken  at  the  same  select  times  as  in 
Figure  4.1.2. 


The  predicted  contours  of  constant  turbulent  fluctuations  are  presented  in  Figure  4.1.1  as  a  function 
of  time  and  altitude.  The  boundary-layer  thickness  grows  during  the  day  and  shrinks  during  the  evening 
and  morning  hours  as  expected.  The  turbulent  kinetic  energy  reaches  a  maximum  of  4.8%  of  the  geostrophic 
mean  kinetic  energy  in  the  midafternoon  at  an  altitude  of  approximately  500  m.  As  the  sun  sets,  the 
turbulence  begins  to  decay  until  in  the  early  morning  hours  the  maximum  kinetic  energy  is  K  0.25*  of  the 
mean  kinetic  energy.  The  biggest  difference  between  the  results  presented  in  Figure  4.1.1  and  those 
presented  by  Lewellen,  et  al.  (1974)  occurs  during  the  early  morning  hours.  The  ^ull  equations  with  a 
dynamic  scale  predict  a  slower  decay  of  q  and  consequently  considerably  larger  q  in  the  altitude  range 
from  100  m  to  1  km  in  the  morning  hours.  But  even  this  larger  value  is  still  quite  small  compared  with 
afternoon  values.  Both  model  representations  predict  such  features  as  the  temperature  inversion  and  local 
wind  jet  observed  to  occur  nocturnally  at  low  levels.  Figures  4.1.2  and  4.1.3  show  the  vertical  profile 
of  the  wind  components  at  selected  times  for  the  lowest  1  km  altitude.  In  the  plots  given  in  Sections  3.4 
and  3.5  for  quasi-steady  conditions,  it  was  possible  to  present  results  in  terms  of  a  dimensionless 
height,  zf/Ug.  However,  in  the  present  case,  since  the  dynamics  of  the  diurnal  variation  are  important  in 
determining  the  height  of  the  layer  at  any  given  time,  this  nondimensional ization  is  no  longer  valid.  For 
Figures  4.1.2  and  4.1.3  the  geostrophic  wind  was  held  constant  at  Ug  =  10  m/sec,  the  upper  level 
temperature  gradient  was  held  constant  at  (*a/*z)|»«  3°C/km,  the  surface  roughness  constant  at 
z0  *  0.01  m,  and  the  surface  heat  flux  was  allowed  to  follow  an  experimentally  observed  distribution 
(Myngaard,  1973).  Between  sunset  and  sunrise,  the  heat  flux  is  a  negative  constant  (-0.025°C  m/sec)  and 
between  approximately  B:0Q  am  and  8:00  pm,  the  heat  flux  is  positive  and  nearly  follows  a  sine  wave  with  a 
maximum  of  0.25°C  m/sec.  Over  the  24-hour  cycle,  there  is  a  net  positive  heat  flux  to  the  atmosphere 
which  must  be  removed  (in  simulation  of  long  wave  radiation)  in  order  to  achieve  a  periodic  solution  with 

3o  long-term  heating  of  the  boundary  layer.  In  the  present  run,  this  was  accomplished  by  assuming  Qr 
istributed  uniformly  over  the  night  hours  with  the  following  assumed  vertical  distribution: 


Qr  *  -  0.703  xlO"1*  cos  (nz/6000m)(°C/sec)  z  <  3000m 


Qr  «  0  z  >  3000  m 


Figures  4.1.2  and  4.1.3  show  the  wind  shear  to  be  a  maximum  in  the  early  morning  hours  when  the 
turbulence,  as  shown  in  Figure  4.1.1,  is  a  minimum.  There  is  a  low-level  nocturnal  jet  which  exceeds  the 
geostrophic  wind  by  14%  at  a  little  above  100  m  altitude.  The  wind  at  this  time  drops  from  11.4  m/sec  to 
3  m/sec  at  10  m  altitude  at  the  same  time  as  the  principal  wind  direction  is  shifting  by  40  degrees. 

Surface  shear  stress  is  plotted  as  function  of  the  stability  variable  *u*/fl  in  Figure  4.1.4.  Data 
points  for  Ro  -  10  ,  as  taken  from  Tennekes'  (1973)  summary,  are  Included  on  the  plot.  The  model 
predictions  show  a  hysteresis  loop  with  the  surface  shear  stress  significantly  larger  when  the  surface 
heat  flux  is  decreasing  rather  than  increasing.  The  data  show  considerable  scatter  but  tend  to  be  biased 
toward  the  upper  bounds.  The  factor  of  2  difference  in  u*  at  neutral  conditions  demonstrates  the 
influence  of  unsteady  dynamics  on  the  atmospheric  boundary  layer. 

A  boundary  layer  data  set  collected  at  Hangara,  Australia  by  Clarke,  et  al.  (1971)  has  been  used  by  a 
number  of  investigators  to  test  their  model  simulations  of  diurnal  variations  in  the  boundary  layer.  The 
data  set  consists  of  the  time  and  height  variation  of  the  two  mean  velocity  components,  the  mean  virtual 
potential  temperature,  and  the  mean  water  vapor  mixing  ratio.  Profiles  of  the  turbulence  quantities  are 
not  Included.  Figure  4.1.5  shows  the  comparison  between  observations  and  computed  results  for  the  diurnal 
variations  of  the  mixed  layer  height  as  given  by  Yamada  and  Mellor  (1975).  The  qualitative  comparison  is 
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Figure  4.1.4  -  Surface  shear  velocity 
vs.  stability  for  diurnal  planetary  boundary 
layers.  Model  calculation  from  Lewellen,  et 
al.  (1974);  data  from  Tennekes  (1973). 
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Figure  4.1.5  -  Time  evolution  of  the  computed  and  observed  mixed  layer  height  as  given  by 
Mellor  and  Yamada  (1975). 


good,  with  the  strong  afternoon  turbulence  responsible  for  the  observed  high  penetration  of  the  mixed 
layer  at  these  times.  Mellor  and  Yamada  report  that  the  approximate  300  m  discrepancy  between  the 
observed  and  computed  boundary-layer  depth  could  probably  be  attributed  to  the  uncertainties  in  the 
observed  values  of  the  mean  vertical  wind  which  served  as  input  to  the  model.  They  concluded  that  fairly 
accurate  data  on  both  the  thermal  wind  and  the  mean  vertical  wind  are  necessary  to  obtain  realistic 
simulations  of  the  mean  variables. 

Andre,  et  al.  (1978)  also  simulated  the  Wangara  data  using  their  third-order  closure  model. 
Comparison  between  the  computed  and  observed  horizontal  wind  profile  during  night  33-34  is  shown  in  Figure 
4.1.6.  The  most  important  feature  is  the  development  of  the  low-level  nocturnal  jet  whose  magnitude 
increases  from  9  m/sec  at  2100  to  12.5  m/sec  at  0300  and  then  decreases.  The  behavior  is  qualitatively 
reproduced  by  the  model  except  the  low-level  jet  is  approximately  100  m  too  low.  They  attribute  this 
discrepancy  to  uncertainty  in  the  aerodynamic  roughness,  zo.  since  the  height  of  the  low-level  jet  is 
dependent  on  the  turbulence  produced  at  the  surface.  They  also  show  that  radiative  heat  transfer  becomes 
more  Important  than  turbulent  transport  above  approximately  100  m  during  nocturnal  hours.  This  will  be 
discussed  in  Section  4.5. 
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Figure  4.1.6  -  Comparison  of  model  calculations 
by  Andre,  et  al.  (1978)  with  Wangara  data. 

Computed  (left)  and  observed  (right)  profiles  of 

mean  horizontal  wind  during  Night  33-34: _ , 

t  =  18  h;  - ,  t  »  21  h;  - ,  t  »  00  h; 

•••*,  t  *  03  h;  —  —  t  =  6  h; _ , 

t  =  08  h. 

Bodin  (1980)  claims  that  fairly  simple  turbulent  closure  models  such  as  his  turbulent  kinetic  energy 
model  and  Long  and  Shaffer's  (1975)  eddy  viscosity  model  do  as  good  a  job  at  simulating  the  Wangara  data 
as  do  the  more  complicated  models.  This  probably  is  more  a  reflection  of  the  importance  of  such  diverse 
phenomena  as  advection,  thermal  wind,  and  radiation  in  this  data  set,  than  it  is  a  verification  of  the 
simple  turbulent  models. 

4.2  Influence  of  Barocl inicity 

The  approximation  of  horizontal,  spatial  homogeneity  is  seldom  completely  realized  in  the 
atmosphere.  One  effect  which  is  often  encountered  is  that  of  horizontal  temperature  gradients. 
Fortunately,  this  effect  can  be  incorporated  within  a  one-dimensional  model  by  permitting  the  pressure 
gradient,  represented  in  Eqs.  (3.2.1  and  3.2.2)  in  terms  of  geostrophic  velocities,  to  vary  with  altitude. 
This  "thermal  wind"  variation  in  the  geostrophic  wind  can  be  obtained  directly  from  the  geostrophic 
balance  by  differentiating  and  using  the  equation  of  state 
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with  x  in  the  direction  of  Ug  at  the  top  of  the  boundary  layer.  Temperature  advection  could  also  be 
included  in  the  energy  equation,  through  a  source  term  of  the  form 
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but  this  is  a  less  important  term.  Neglecting  temperature  advection,  this  introduces  two  parameters  into 
the  governing  equations  to  specify  the  magnitude  and  direction  of  the  geostrophic  wind  shear  (Arya  and 
Wyngaard,  '975).  These  will  be  taken  herein  as 
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Arya  and  Wyngaard  (1975)  show  that  when  M  and  do  are  steady  and  the  boundary  layer  is  unstable,  then 
most  of  the  imposed  thermal  wind  shear  occurs  at  the  top  of  the  boundary  layer.  Even  for  only  slightly 
unstable  conditions  the  mean  velocity  components  remain  essentially  uniform  across  the  “mixed"  layer  with 
the  wind  shear  occurring  across  or  above  the  capping  inversion  layer.  Under  these  conditions,  the 
principal  influence  of  adding  the  thermal  wind  is  to  increase  the  wind  shear  across  the  inversion  layer. 
This  should  tend  to  enhance  the  potential  for  shear-wave  instabilities  above  the  inversion  and  a  general 
thickening  of  this  layer. 

Conditions  which  lead  to  significant  values  of  M  in  the  boundary  layer  are  generally  not  steady. 
Lewellen  and  Williamson  (1976)  show  that  whenever  the  time  scale  for  the  development  of  barocl inicity  is 
of  the  order  of  a  few  hours,  as  it  is  likely  to  be  in  nature,  then  the  time  history  of  the  barocl inicity 
affects  the  boundary -layer  response.  Figures  4.2.1  and  4.2.2  show  results  for  a  constant  value  of 
M  (1.25)  and  periodic  change  in  i0  between  0  and  2*  over  a  period  of  12  hours.  The  magnitude  of  M 
corresponds  to  a  horizontal  temperature  change  of  3“C  in  100  km,  a  large  but  not  uncommon  value  in  nature. 
The  time  variation  corresponds  to  a  relatively  fast  moving  temperature  anomaly  with  the  flow  permitted  to 
become  periodic  by  allowing  three  full  periods  of  variation  before  the  results  shown  were  recorded.  The 
geostrophic  wind  was  set  at  10  m/sec,  Ro  =  10“  and  no  surface  heat  flux,  Ri(io)  =  0.  An  upper  level  lapse 
rate  of  3°C/km  held  the  boundary-layer  thickness  to  «  1280  m  altitude. 

U/Ug  CONTOURS 


Figure  4.2.1  -  Isolines  of  U  for  a  thermal  wind  varying  periodically  in  direction 
[Ro  *  10“,  au/az|»  =  3°C/km,  M  *  1.2b,  Ri^yj  *  UJ.  Lewellen  and  Williamson  (1976). 
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Figure  4.2.2  -  Isolines  of  V  for  a  thermal  wind  varying  periodically  in  direction 
[Ho  -  10“ ,  aa/vz  |  •  ■  J°C/km,  M  ■  1.25,  Ri^yj  •  0].  Lewellen  and  Williamson  (1976). 


From  Figure  4.2.1  it  is  seen  that  the  maximum  U  velocity,  13.5  m/sec,  occurs  approximately  3  hours 
after  the  thermal  wind  is  directed  in  the  direction  of  the  geostrophic  wind  (80  =  0).  however,  3  hours 
after  the  thermal  wind  is  next  directed  in  the  same  direction  (30  -  2x),  U  is  near  a  minimum.  Actually, 
from  Eq.  (3.2.1)  and  Eq.  (3.2.2),  we  see  that  the  U  velocity  responds  to  variations  in  Vg  and,  since 
3Vg/3z  follows  sin  3t,  this  contributes  a  positive  pressure  gradient  between  0  <  a0  <  *  and  a  negative 
pressure  gradient  between  *  <  80  <  2*.  Thus  it  is  consistent  for  a  maximum  to  occur  in  11  when  n0  passes 
through  ”  as  it  is  increased,  and  for  a  minimum  to  occur  in  U  when  30  passes  through  ■  as  it  is  decreased. 
This  particular  time  variation  of  driving  conditions  is  not  expected  to  occur  long  enough  in  nature  to 
make  the  real  flow  periodic,  but  it  does  demonstrate  how  important  the  time  history  of  the  flow  is  in 
determining  the  response  of  the  flow  to  instantaneous  values  of  the  driving  conditions. 

The  maximum  value  of  V,  7  m/sec,  occurs  near  the  ten-hour  mark  after  the  combination  of  a  relative 
minimum  in  U  and  a  maximum  in  Ug  has  produced  a  favorable  pressure  gradient  [see  Eq.  (3.2.2JJ  for  a  period 
of  time.  There  is  a  region  of  slightly  negative  V  (approximately  -1  m/sec)  between  1000  m  and  1230  in 
altitude  which  does  not  show  in  Figure  4.2.2.  The  turbulence  isolines  (Lewellen  and  Williamson,  1976) 
show  no  strong  variation  over  the  period. 

The  maximum  wind  shear  at  or  above  10  m  altitude  computed  by  taking  the  square  root  of  the  sun  of  the 
squares  of  3U/3z  and  <*V/«*z  gives  a  value  of  0.14  sec'1,  somewhat  less  than  the  maximum  value  (0.18  sec'1) 
that  occurs  in  the  typical  diurnal  day  calculations  of  the  Section  4.1.  However,  in  this  case  the  shear 
occurs  over  a  somewhat  wider  altitude  range.  In  the  earlier  figures  the  strong  local  stability  caused  by 
surface  cooling  was  the  dominant  mechanism  while,  in  the  present  case,  there  is  a  neutral,  vertical 
temperature  variation  in  the  region  of  maximum  shear.  This  exemplifies  two  different  ways  in  which  the 
wind  shear  may  be  increased:  either  by  increasing  the  horizontal  wind  gradients  to  increase  the  driving 
thermal  winds,  or  by  applying  a  vertical  temperature  gradient  which  increases  the  stability  of  the  flow 
under  the  same  driving  pressure  gradients. 

In  any  attempt  to  simulate  conditions  in  the  baroclinic  atmospheric  boundary  layer  at  a  given  time 
and  place,  it  appears  desirable  to  trace  the  development  of  the  boundary  layer  for  some  time  prior  to  the 
desired  time. 


4.3  The  Marine  Boundary  Layer 

The  marine  boundary  layer,  recently  reviewed  by  Lemone  (1980),  differs  in  several  respects  from  that 
over  land.  The  sea  surface  tends  to  be  quite  uniform  with  relatively  slow  changes  in  temperature. 
Although  there  are  exceptional  areas  where  regions  of  up welling,  special  currents,  or  coastal  influence 
may  permit  horizontal  temperature  gradients  of  several  degrees  over  a  few  kilometers,  temperature 
differences  over  spatial  dimensions  of  the  order  of  tens  of  kilometers  and  over  diurnal  time  scales  are 
typically  less  than  a  degree.  The  sea  is  a  fluid  lower  boundary  which  interacts  with  the  atmosphere  to 
determine  in  a  coupled  manner  the  effective  aerodynamic  roughness  of  the  surface.  Also  the  marine 
boundary -layer  air  is  typically  quite  moist.  In  fact  the  surface  buoyancy  flux  is  often  determined  as 
much  or  more  by  the  moisture  flux  as  it  is  by  the  temperature  flux. 

For  the  marine  boundary  layer,  it  is  desirable  to  replace  the  potential  temperature  by  the  virtual 
potential  temperature  defined  as 


uy  «  e  +  0.61  T0  h 


(4.3.1) 


As  long  as  there  is  no  condensation  or  evaporation  of  water  within  the  boundary  this  replacement  permits 
all  the  previous  equations  to  be  applied  to  the  marine  boundary  layer.  Change  of  phase  will  be  considered 
in  Section  4.4. 

The  air-sea  surface  interaction  could  rightfully  be  the  subject  of  another  monograph  (see,  e.g.. 
Coantic,  1978).  Many  details  of  this  interaction,  particularly  under  conditions  of  significant  sea  spray, 
have  not  been  quantitatively  analyzed.  The  simplest  relation  with  some,  at  least  approximate,  general 
validity  is  to  use  Charnock's  (1955)  Froude  number  scaling  to  relate  the  effective  surface  roughness  to 
the  surface  shear  stress.  When  Wu’s  (1969)  data  correlation  is  used  to  evaluate  the  constant,  this  gives. 


z0  =  0.016  Ui/g 


(4.3.2) 


The  type  profiles  which  typically  occur  in  the  marine  boundary  layer  are  shown  in  Figure  4.3.1  from 
Lemone  (1980).  A  slightly  unstable  mixed  layer  capped  by  a  weak  inversion  is  common.  Comparison  between 
the  heat  flux  and  moisture  flux  measured  by  aircraft  on  one  of  these  days  and  the  results  of  calculations 
by  the  model  described  by  Eqs.  (2.3.18  -  2.3. 23)  are  given  in  Figure  4.3.2  as  given  by 
Lewellen,  et  al.  (1977).  These  comparisons  were  used  to  show  that  the  model  is  roughly  consistent  with 
the  data.  A  careful  assessment  of  the  accuracy  of  the  model  would  require  a  more  detailed  specification 
of  the  experimentally  observed  boundary  layer  as  it  evolves  in  time. 

The  most  sophisticated  attempt  at  calculating  the  marine  boundary  layer  is  that  by  Sommeria  (1976) 
and  Sommeria  and  Lemone  (1978),  using  a  three-dimensional  simulation  of  the  turbulent  processes  in  a 
manner  similar  to  the  calculations  of  Deardorff  (1972,  1974).  The  biggest  difference  between  their  model 
results  and  that  of  the  one-dimensional,  second-order  closure  result  of  Lewellen,  et  al.  (1977),  appeared 
to  be  a  significant  difference  in  liquid  water  content  condensed  in  the  clouds.  The  role  of  moisture 
change  of  phase  is  quite  important  in  the  marine  boundary  layer.  This  is  discussed  in  Section  4.4. 
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Figure  4.3.1  -  Profile  of  potential  temperature,  t>;  virtual  potential 

and  specific  tumidity,  q;  for  4  days  in  GATE,  1974  (Day  218  =  6  Aug.; 
Day  253  =  10  Sept.;  Day  248  =  15  Sept.  (LeMone,  1980). 


temperature,  ovi 
Day  243  =  31  Aug.; 


IS  Sept  GATE  CASE 
Reference  time  •  12  .00  hr 


Figure  4.3.2  -  Virtual  potential  temperature  and  humidity  fluxes  as  calculated  by  a 
second-order  closure  model  (Lewellen,  et  al.,  1977),  and  compared  with  the  GATE 
measurements  (LeMone,  1980). 


4.4  Moisture  Change  of  Phase 

The  latent  energy  transported  through  the  boundary  layer  in  the  form  of  water  vapor  is  ultimately 
responsible  for  fueling  most  of  the  large  scale  atmospheric  motions.  On  a  global  average  basis,  the 
turbulent  transport  of  water  vapor  through  the  atmospheric  boundary  layer  must  balance  the  average 
precipitation  of  liquid  water,  about  1  m  per  year  (Coantic,  1978).  Even  if  only  a  small  fraction  of  the 
energy  released  by  condensation  of  the  water  vapor  occurs  directly  in  the  boundary  layer,  this  is 
1  t0  'n^era£t.  strongly  with  the  turbulent  structure  of  the  boundary  layer  when  it  does  occur. 
The  other  reason  for  being  interested  in  the  interaction  of  turbulence  with  moisture  change  of  phase  in 
the  boundary  layer  is  that  the  phenomena  of  fog  and  low-level  clouds  represent  a  primary  goal  of  many 
models. 
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In  the  presence  of  moisture  change  of  phase,  the  potential  temperature  is  not  always  the  most 
convenient  variable  to  use  in  the  energy  equation,  since  it  is  not  a  conserved  quantity.  One  choice  for  a 
conserved  quantity  is  the  liquid  potential  teinperature  (Betts,  1973,  1974). 


T 


To 


(4.4.1) 


An  alternative  is  the  moist  static  energy. 


To 


(4.4.2) 


The  difference  between  ©s  and  o«,  is  the  latent  energy  of  the  total  water  content.  As  long  as  the  total 
water  content  is  also  a  conserved  quantity,  the  two  are  equivalent.  In  the  presence  of  precipitation,  ©s 
is  conserved  but  o*  is  not. 

The  buoyancy  term  in  the  momentum  equation  is  proportional  to  the  virtual  potential  temperature 
defined  as 


©v  *  T  -  T0  +  ~  z  *  (0.61Hv-H*)To 
CP 

Thus 


(4.4.3) 


(4.4.4) 


so  that  even  with  both  os  and  H  conserved,  a  change  of  phase  of  part  of  H  from  Hy  to  Ht  leads  to  a  source 
of  ©y. 

Because  of  the  extensive  role  ©y  plays  in  the  turbulent  flux  equations  [Eqs.  (2. 3. Id  -  2.3.21)  with 
©v  replacing  0],  Oliver,  et  al .  (1978)  chose  to  stick  with  ©y  as  a  primary  variable  and  introduce  the 
necessary  source  terms.  The  energy  equation  is  then  replaced  by 


D<»y  3  _  .  g 

Ot  =  '  3Xj  {uJB)  '  Ks  grad  +  K9  Cp  w  6iJ 


(4.4. b) 


where  K$  and  Kg  are  coefficients  generated  from  the  vapor  saturation  function 


(1  ♦  1.61HS  dT) 

I  Cl  ♦  (Ry/Cp)  Hs*tJ 


1  ♦  1.61  Hs  ^ 
S  Rm 


(1  ♦  1.61  Hsar) 

(1  *  (Ry/Cp)  H$  dfj 


(4.4.7) 


The  upper  values  of  Ks.Kg  apply  for  unsaturated  conditions  and  the  lower  values  for  saturation  conditions. 
Ry  and  Rm  are  gas  constants  for  water  vapor  and  air,  respectively.  The  parameter  dj  is  the  logarithmic 
derivative  of  the  saturation  mixing  ratio  with  respect  to  temperature.  If  the  saturation  formation  is 
given  by  the  Clapeyron  relationship 
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where  the  subscript  “denotes  reference  values,  then  By  is  given  by 


3T  =  sL 

T  RuT 


(4.4.8) 


(4.4.9) 


Additional  terms  also  need  to  be  added  to  the  heat  flux,  temperature  variance,  and  temperature-humidity 
covariance  equations.  These  are: 


..  +  Kg  U1w 


Dt  =  —  +  2K9^W°v 
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These  additional  terms  play  the  role  of  reducing  the  effective  vertical  temperature  gradient  in  the 
production  terms  in  these  equations,  i.e.,  the  effective  vertical  temperature  gradient  appearing  in  each 
of  these  equations  is  [aey/az  -  Kg(g/cn)].  The  presence  of  intense  turbulence  will  force  the  ambient 
environment  towards  a  lapse  rate  of  (Kg  -  l)g/cp,  which  is  the  moist  adiabatic  lapse  rate  under  saturated 
conditions,  or  the  dry  adiabatic  lapse  rate  under  unsaturated  conditions.  Just  as  the  dry  adiabatic  lapse 
rate  (~  -9.8  °C/km)  marks  the  boundary  between  stable  and  unstable  considerations  in  a  dry  environment, 
the  moist  adiabatic  lapse  rate  -3.5  oc/km)  marks  the  stability  boundary  under  saturated  conditions.  If 
the  ambient  virtual  potential  temperature  gradient  falls  between  these  two  stability  boundaries,  then  a 
conditional  instability  exists.  As  long  as  the  flow  remains  unsaturated,  the  flow  will  remain  stable,  but 
if  there  is  a  disturbance  of  sufficient  amplitude  to  allow  saturation  to  be  reached  at  some  point,  then 
this  disturbance  will  become  unstable. 

Just  as  moisture  change  of  phase  has  a  strong  input  on  buoyantly  generated  turbulence,  turbulence  may 
be  expected  to  have  a  strong  interaction  with  the  microphysics  of  phase  change.  A  discourse  on  cloud 
microphysics  is  beyond  the  scope  of  the  present  monograph  (e.g.,  see  Pruppacher  and  Klett,  1978),  although 
I  believe  it  is  an  area  where  turbulence  modeling  could  fruitfully  be  applied.  The  simplest  approximation 
is  to  assume  equilibrium  condensation  so  that  the  vapor  mixing  ratio  Hv  and  the  liquid  mixing  ratio  Ht  are 
given  by 


Hv  +  hv  =  H  +  h  +  (Hs  +  hs-H-h)fcRH  +  h-Hs-hs) 


Ht  +  ht  =  (H  +  h-H,+hs)#(H  +  h-Hs-h,) 


(4.4.13) 
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when#- is  the  Heaviside  function  and  Hs(p,T)  the  saturation  mixing  ratio  given  by  Eq.  (4.4.8). 
It  is  fairly  apparent  that  the  ensemble  average  liquid  water  content  given  by 


is  not  equal  to  the  simple  form 


*  (H  +  h  -  Hs  +  hs)£f-(H  +  h  -  Hs  -  hs) 


(H-HS)#(H-HS) 
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because  liquid  water  content  is  physically  a  positive  definite  quantity.  The  exact  specification  of  the 
average  liquid  water  content  depends  on  the  joint  probability  distribution  factor,  f(H,Hs).  Sommeria  and 
Oeardorff  (1977)  and  Hellor  (1977)  assumed  a  joint  Gaussian  distribution  about  the  mean  for  the  2 
conserved  quantites  ®s  and  H  to  obtain  a  result  which  may  be  parameterized  as 


Qz  <  -  1.6 


Ht  *  v  02  »  (qz  +  1.6)V6.4  -1.6  <Q2  <1.6 


Qz  »  1.6 
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where 


°2  *  Clh-hs)2]  1/2 


and 


Qz  =  (H-Hs)/o2 


Oliver,  et  al.  (1978)  adopted  a  somewhat  more  approximate  procedure.  They  approximated  the  average  liquid 
water  content  as 


Hi 


H  -  Hs 

\  (H  -  Hs  +  a) 


a  <  |H-HS| 
a  >  |H  -  Hs| 


(4.4.17) 


where  a  is  a  characteristic  measure  of  the  size  of  the  fluctuations  which  may  be  partially  satutaced  and 
partially  unsaturated: 


a 


(h2) 


1/2 


(4.4.18) 


The  functions  tU/<>2  and  Hi/a  are  quite  similar  functions  of  Q2  and  [H  -  Hs)/“,  but  a  >  02  so  that  the 
latter  approach  leads  to  a  broader  transition  region  and  more  liquid  water  content  over  most  of  the 
transition  layer. 

Water  droplets  are  diffused  by  turbulence  in  the  same  manner  as  a  gaseous  species  as  long  as  the 
particle  is  sufficiently  small  that  it  can  faithfully  follow  all  of  the  turbulent  eddies.  This  will  be 
true  whenever  the  relaxation  time  of  the  particle  tv,  the  time  required  by  a  particle  to  reduce  its 
velocity  relative  to  the  gas  by  e-1  when  acted  upon  by  the  fluid  drag  only,  is  much  less  than  the  eddy 
time  scale  A/q,  Table  4.4.1  gives  some  values  of  ty  for  different  size  particles  with  a  density  of 
1  gm/cm*.  Since  A  grows  proportional  to  z  near  the  surface,  it  is  possible  to  estimate  the  height  at 
which  qty/A  =  1  for  a  given  value  of  q  (1  m/sec).  At  large  distances  above  this  height,  there  should  be 
no  difference  between  the  turbulent  diffusion  of  a  particle  and  the  turbulent  diffusion  of  a  gaseous 
species.  For  cloud  droplets  with  diameters  less  than  approximately  10  microns,  the  table  indicates  that 
the  only  region  where  the  particulate  nature  will  generally  make  any  difference  to  diffusion  is  within  the 
surface  sublayer.  When  the  droplet  stopping  distance  is  of  the  order  of  the  turbulent  microscale,  or 
larger,  efficient  deposition  of  fog  droplets  on  the  surface  may  be  expected.  It  can  also  be  seen  from 
Table  4.4.1  that  the  gravitational  settling  contribution  to  deposition  velocity  will  only  be  important  for 
droplets  of  the  order  of  10  wn  or  larger. 


Table  4.4.1.  Aerosol  Transport  Properties 

(Spherical  particles  with  «  1  gm/cm*  in  air  at 
20°C,  1  ATM;  and  q  *  1  m/sec) 


Particle 

diameter 

2a 

(um) 

Relaxation 

time 

Tyt 

(sec) 

Terminal 

velocity 

Tv9 

(cm/ sec) 

Height  at 
which 

q*iy/A  •  1 

(m) 

Stopping  distance 
microscale 

q2vv/2v 

0.1 

8.8  *  10- 8 

8.6  x  10-s 

2.9  x  10-3 

0.5 

1.0  *  10-8 

1.0  x  10-* 

3.3  x  10-4 

1.0 

3.6  *  10-8 

3.5  x  10-3 

0.12 

5 

8.0  «  10-s 

7.8  x  10-* 

12.3  x  10-s 

2.7 

10 

3.2  x  10-8 

3.1  x  10- 1 

5  x  10-" 

10. b 

50 

7.7  x  10-3 

7.58 

0.012 

2.6  x  104 

100 

3.1  x  10-* 

30.3 

0.048 

1.0  x  101 

t  Friedlander  (1977) 
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4.5  Thermal  Radiation 

The  principal  mode  of  transport  of  energy  through  the  atmosphere  is  via  thermal  radiation,  but  the 
relative  transparency  of  the  general  atmosphere  allows  this  energy  to  typically  pass  through  a  layer  of 
the  order  of  1  km  thick  without  a  strong  interaction.  Of  course,  the  thermal  radiation  balance  always 
plays  an  important  role  in  controlling  the  surface  temperature  and  thus,  at  least  indirectly,  plays  a 
deciding  role  in  determing  the  stability  of  the  boundary  layers  as  discussed  in  Chapter  III.  However, 
under  some  conditions,  particularly  those  involving  fog  or  low-level  clouds,  thermal  radiation  absorption 
depths  may  be  greatly  reduced  and  permit  a  strong  interaction  between  turbulent  transport  and  thermal 
radiation  within  the  atmospheric  boundary  layer.  Also  under  nocturnal,  clear  sky  conditions,  the 
radiation  flux  divergence  may  be  stronger  than  the  divergence  of  the  turbulent  heat  flux.  Thus,  a  general 
model  of  the  atmospheric  boundary  layer  needs  to  at  least  include  thermal  radiation  in  the  mean  energy 
equation.  The  extent  to  which  radiative  terms  need  to  be  included  in  the  temperature  correlation 
equations  remains  unresolved.  Andre,  et  al.  (1978)  argue  that  a  radiative  dissipation  term  in  the 
temperature  variance  equation  is  important  for  simulating  the  nocturnal  boundary  layer  under  even  clear 
conditions  as  long  as  there  is  sufficient  water  vapor  present. 

In  fornulating  a  description  of  radiant  heat  transfer,  it  is  convenient  to  separate  the  direct  solar 
(short  wavelength)  from  the  terrestrial  (long  wavelength)  fields.  Each  of  these  fields  may  be  described 
in  terms  of  the  one-dimensional  (z)  upward  and  downward  frequency  and  angular  averaged  radiation 
intensities  F+,  F".  These  intensities  are  governed  by  the  upward  and  downward  transport  equations  in 
integral  form  (Goody,  1964): 


F+(z)  *  [F+(z0)  -  l(z0)]  t(z0,z)  +  i(z)  -  j  TU'.zl^dz1  (4.5.1) 

*o 

z 

F*(z)  *  [F-(zh)  -  i(th)3  *(zh.z)  +  t(z)  -  f  t(z‘ ,z)  dz‘  (4.5.8) 

Zh 

The  boundary  values  at  the  surface  z  ■  z„  and  boundary  layer  “top"  z  -  are  F +(z0)  and  F*(zh), 
respectively.  The  radiative  cooling  rate,  Qr  Is  given  in  terms  of  F+,  F*  for  either  radiation  field  as 


(Jr  -  (F+  -  F*)  (4.5.3) 

oZ 


and  the  total  radiant  heat  flux  is  the  sum  of  the  direct  solar  heat  flux  and  the  terrestrial  heat  flux. 

The  absorption,  scattering,  and  emissive  properties  of  the  atmosphere  must  be  embedded  in  the  source 
function  i(z)  and  the  transmission  function  tjzi.z^,)  between  two  levels  Zj.Zy.  A  simple  analysis 
neglecting  scattering  of  direct  solar  radiation  was  given  by  Oliver,  et  al.  (1978).  Their  analysis 
demonstrates  that  the  influence  of  liquid  water  on  the  transmission  function  x  is  quite  different  for 
solar  radiation  than  It  Is  for  terrestrial  radiation.  For  solar  radiation  the  transmission  function  x  is 
modeled  after  Manabe  and  Strickler  (1964)  for  water  vapor  and  CO*.  To  account  for  liquid  water,  a 
transmission  xfiq(zj,z2)  was  included  as  a  factor  in  the  total  transmission  x  of  the  form 


xfiqUi.zJ  «  exp(-  iy  usm)  (4.5.4) 


where  m  Is  the  absorbing  mass 


(4.5.5) 


6  is  the  dlffuslvity  factor,  and  us  Is  the  specific  absorption  coefficient  for  direct  solar  radiation  by 
water  droplets.  Feigel'son  (1964)  notes  the  solar  absorption  is  not  sensitive  to  droplet  size  for 
wavelengths  <  2  un.  Integrating  the  spectral  liquid  water  absorption  coefficient  tabulated  by  Feigel’son 
over  the  direct  solar  spectrun,  gives 


16  cmzg*1 


For  the  terrestrial  field,  the  transmission  function  for  water  vapor  and  CO?  follows  that  of 
Feigel'son  (1970)  with  the  liquid  water  transmission  included  as  a  factor  in  the  tranmission  i  of  the  form 


exp(-  iiu  m) 


(4.5.6) 


Here  is  the  averaged  liquid  water  absorption  coefficient  for  terrestrial  radiation  which  depends 
strongly  upon  the  droplet  size  distribution  ranging  in  value  from  200  to  1700  cm^g-1  (Feigel'son,  1970) 
for  spectral  width  average  drop  size  radii  ranging  from  4.5  to  7  wn.  The  droplet  size  effect  on  radiative 
absorption  may  be  maintained  as  a  parameter  through  uT.  Detailed  inclusion  of  variation  with  droplet  size 
requires  prediction  of  the  droplet  number  density  within  a  cloud  or  fog  as  well  as  the  liquid  water 
content.  Experimental  evidence  (Mack,  et  al.  1974)  does  not  validate  the  simplifying  assumption  of 
uniform  droplet  density  in  space.  Their  observations  indicate  an  order  of  magnitude  variation  in  droplet 
number  density  through  the  boundary  layer  while  mean  droplet  size  varies  by  about  30*.  A  characteristic 
value  of  J  for  fog  and  stratus  (mean  drop  radius  ~  6  ym)  is  600  cm^g-1. 

The  significant  difference  between  the  absorption  coefficients  us  and  J  is  noteworthy.  The 
absorption  coefficents  may  be  represented  in  terms  of  absorption  lengths  as 


x 


T 

rad 


*rad 


us  u  <H*> 


(4.5.7) 


which  are  exhibited  in  Table  4.5.1  for  various  (spatial)  average  liquid  water  contents  (Hx).  These  values 
indicate  that  water  is  radiatively  more  active  with  the  terrestrial  spectrum  than  the  direct  solar 
spectrum  for  drops  in  the  4-10  wn  range  typical  of  fog  and  cloud.  This  sharp  difference  in  radiative 
emission/absorption  in  the  terrestrial  and  solar  spectrum  influences  strongly  the  structure  of  cloud/fog 
cooling  and  the  mechanism  by  which  direct  solar  radiation  evaporates  away  clouds  or  fog.  These 
implications  are  discussed  in  the  following  sections. 

Appropriate  boundary  conditions  for  terrestrial  radiative  transport  are  F*  =  0  at  the  upper  boundary, 
and  F+  *  B0  (B0  is  the  blackbody  flux)  at  the  surface.  For  direct  solar  radiation,  F"  is  specified  as 
the  incident  downward  solar  flux  at  the  upper  boundary  and  F+  =  rF‘  at  the  surface,  where  r  is  the  surface 
reflectance.  If  the  upper  boundary  for  radiative  transport  is  set  at  the  atmospheric  scale  height,  the 
average  liquid  water  content  and  average  water  vapor  content  of  the  atmosphere  contained  within  the 
radiative  boundary  layer  above  the  top  of  the  atmospheric  boundary  layer  must  be  specified. 


Table  4.5.1. 


Radiative  absorption  depths  for  fogs  and  stratus. 


,s,T 

Tad 


s,T 


Ho  J 


-1 


Average 
liquid  water 
contents 

Hi 

(g  kg-i) 

Direct  solar 
radiation 

*rad 

(m) 

Terrestrial 

radiation 

*rad 

(m) 

Tenuous 

fog 

0.025 

20,000 

600 

0.1 

5,500 

167 

0.2 

2,700 

83 

0.4 

1,350 

41 

Dense 

0.8 

700 

21 

stratus 

1.0 

500 

16 

tops 

4.6  Fog 

Shallow  fogs  in  which  turbulent  transport  dominates  radiative  transport  may  be  modeled  as  a 
straightforward  extension  of  the  surface  layer  analysis  considered  in  Section  3.1.  In  the  absence  of 
radiation  and  precipitation,  the  equations  for  0s  and  H  are  similar  to  one  another;  both  their  mean 
equations  and  the  turbulent  correlation  equations.  Thus,  the  distribution  of  0s  -  0so  and  H  -  H0  as  well 
as  their  turbulent  correlations  should  be  similar  within  the  surface  layer,  and  functions  of  z/L  and  z/z0 
only.  As  essentially  noted  by  Taylor  (1917),  all  fluid  states  in  the  surface  layer  should  map  onto  a 
straight  line  In  the  H  -  0s  plane  which  passes  through  the  surface  (Ho,0so).  and  has  a  slope  determined  by 
the  ratio  of  the  surface  humidity  flux  to  the  heat  flux,  H./Os*. 


J 
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For  saturated  conditions  at  the  surface,  Oliver,  et  al.  (1978)  showed  that  surface  fog  is  only 
possible  for  a  cold  surface. 


tts*  <  0  ,  If  >WoSs  >  ^  (4,6.1) 


or  for  a  warm  surface 


°s,  >  0  »  1f  <  d0  (4.6.2) 


where  do  is  the  thermodynamic  derivative  with  respect  to  0$  of  the  saturation  curve  Hs. 
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(4.6.3) 


The  fog  regimes  defined  by  Eqs.  (4. 6. 1-4. 6. 3)  are  exhibited  in  Figure  4.6.1  as  a  function  of  surface 
temperature  and  the  relative  humidity  and  temperature  of  the  air  at  some  distance  above  the  surface. 

An  illustrative  family  of  fog  distributions  as  given  by  the  similarity  functions  Eqs.  (3.1.9)  and 
(3.1.11)  is  shown  in  Figure  4.6.2.  Of  particular  interest  is  the  structure  of  the  liquid  water  profiles 
in  the  fog  layer  for  stable  (c  >  0)  and  unstable  (c  <  0)  regimes.  Stable  fog  layers  have  maximum  H»  close 
to  the  surface,  and  drop  off  sharply  in  liquid  water  content  up  to  their  maximum  depth.  Unstable  layers 
are  well-mixed  and  do  not  show  a  pronounced  maximum  except  for  very  thin  layers  a/L  «  1.  Both  families 
of  fogs  deepen,  and  their  maximum  liquid  water  contents  increa:e  as  their  flux  ratios  H«/0s*  depart  more 
from  the  critical  value  (So¬ 
under  what  conditions  will  the  surface  fog  layer  be  influenced  by  radiative  as  well  as  turbulent 
transport?  The  surface  turbulent  heat  flux  will  be  determined  by  the  temperature  gradient  in  the  surface 
layer;  but  this  gradient  will  also  determine  the  maximum  liquid  water  developed  in  the  fog  layer  and 
hence,  the  radiative  emission  from  the  layer.  For  a  radiatively  thin  fog  layer  (6/lrad  «  1),  it  is 
possible  to  estimate  the  radiative  to  turbulent  heat  flux  for  the  fog  layer  as 


(4.6.4) 


For  a  nominal  condition  of  a  300  K  surface,  [u0  -  H*/es,]  z  2  *  10-5,  and  uT  =  500  cm-Vgm, 

A  9  *  TdO  meters)  ’T(*eo  surface) 
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Figure  4.6.1  -  Schematic  showing  the  air  relative 
humidity  and  air-sea  temperature  differences  required 
for  pure  advective  fogs. 
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Figure  4.6.2  -  Illustrative  liquid  water  distribution  in  stable  ( v  >  0)  and  unstable 
(v  <  0)  surface  layer  advective  fogs,  Oliver,  et  al.  (1978). 


grad  _  a  &n(6/z0) 

qturb  ~  2000u* 


a  result  which  depends  only  upon  the  strength  of  the  wind  and  the  depth  of  the  fog.  It  takes  a  relatively 
strong  wind  to  make  qrad  negligible  compared  to  qturb.  For  6  •  10  m,  and  z0  *  ID'3  m,  a  surface  shear 
velocity  of  u*  =  0.6  m  s-1  (corresponding  to  a  neutral  wind  at  10  m,  of  14  ms-1)  or  greater  is  required  to 
make  qrad  less  than  10%  of  qturb.  Thus,  radiation  plays  a  significant  role  in  most  surface  fogs  of  depths 
greater  than  a  few  meters. 

It  is  necessary  to  return  to  the  full  turbulent-radiative  description  to  illustrate  the  general  case 
of  fog  formation  with  both  radiant  and  turbulent  transport.  Oliver,  et  al.  (1978)  used  the  model 
described  in  Sections  4.3  through  4.5  to  predict  the  formation  and  structure  of  a  comnonly  observed  class 
of  coastal  fogs  over  water  in  which  convection  and  radiation  play  equally  important  roles.  Fogs  are 
frequently  found  to  occur  in  coastal  regions  subject  to  stong  upwelling  and  nonuniform  surface 
temperatures  where  winds  bring  the  surface  air  from  a  region  of  cold  water  onto  a  region  of  warmer  water 
(warm  surface  fog).  Such  fogs  have  been  reported  and  documented  by  Mack,  et  al.  (1973)  and  Leipper 
(1948).  A  stable,  steady  boundary  layer  running  from  constant  temperature  water  to  a  region  of  linearly 
increasingly  surface  temperature  provides  an  archetypal  example  of  such  advective-radiative  fogs.  As  a 
result,  an  unstable  mixed  layer  will  be  initiated,  and  fog  may  form  in  the  newly  developing  layer. 

The  incident  profiles  for  this  illustration  are  characterized  by  a  flux  ratio  of  *  0.0001b 
under  conditions  where 


- j  *  0.0006  . 

Bv*/crit 


The  full  profiles  are  determined  by  a  wind  velocity  of  2.8  ms-1  and  an  absolute  temperature  difference  of 
+1°C  at  z  »  20  m  elevation  above  the  sea  surface  of  T  =  IS'C.  Thus,  the  Incident  flow  is  in  the  cold 
surface  fog-free  regime.  Direct  solar  radiation  was  suppressed  (simulating  nighttime  conditions).  The 
surface  humidity  was  maintained  at  its  saturation  value  at  the  surface  temperature.  Surface  temperature 
linearly  increases  from  the  origin  by  0.5°C  in  3  km.  The  response  of  this  layer  is  shown  in  Figure  4.6.3. 
The  combined  action  of  increasing  surface  temperature  and  radiative  cooling  generates  a  fog  layer  which 
grows  as  the  air  is  convected  over  the  wanner  water.  Near  the  origin  of  the  change  in  surface  temperature 
(x  *  100  m),  a  thin  fog  layer  of  about  18  m  depth  exists  with  maximum  liquid  water  content  just  above  the 
surface.  At  1000  m  the  layer  is  100  m  deep.  At  this  depth,  the  liquid  water  content  has  developed  two 
relative  maxima  -  one  near  the  surface  and  one  near  the  fog  top.  This  is  because  the  100  m  deep  fog  band 
is  no  longer  radiatively  thin,  and  fog-top  radiative  cooling  is  beginnning  to  dominate  over  surface 
radiative  cooling.  The  early  stages  of  a  temperature  inversion  formation  at  the  top  of  the  fog  band  are 
also  in  evidence. 

After  the  change  from  stable  to  unstable  conditions,  and  the  onset  of  radiative  cooling,  the  local 
H*/«v*  shifts  to  a  value  of  0.0008.  Such  a  flux  ratio  places  a  turbulent  dominated  layer  in  the  warm 
surface  fog-free  regime.  Fog,  nevertheless,  exists  in  such  a  layer  because  of  the  radiative  cooling 
present.  A  convenient  measure  of  the  effect  of  radiation  In  a  fog  layer  of  thickness  6  is  an  effective 
radiative  flux  ejad  defined  as 


u.  aCad 
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Direct  comparison  of  ejad  Mi th  es«  provides  an  assessment  of  the  strength  of  radiation  transport  compared 
to  turbulent  transport  in  the  fog  layer.  For  this  warm  surface  advection  fog,  ~  o.b.  The 
role  of  radiation  is  therefore  pivotal  in  its  evolution. 

The  structure  of  radiative  cooling  in  the  fog  is  shown  in  Figure  4.6.4.  In  the  first  1000  m,  maximum 
radiative  cooling  is  concentrated  at  the  surface,  but  then  begins  to  shift  upward  into  the  center  and  top 
of  the  cloud  bank  as  it  becomes  radiatively  thick  to  terrestrial  radiation.  The  turbulence  level  depicted 
in  Figure  4.6.5  shows  a  steady  increase  with  x  as  the  unstable  mixed  layer  grows  with  the  maximum 
turbulence  levels  concentrated  in  the  center  of  the  cloud  bank. 
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Figure  4.6.3  -  Temperature  and 
liquid  water  content  profiles  in  an 
advective-radlatlve  fog  at  x  *  100  m 
(a),  x  =  1000  m  (b),  and  x  *  3000  m 
(c).  Oliver,  et  al.  (1978). 


Figure  4.6.4  -  Radiative  cooling 
distribution  in  advective-radiative 
fog.  Note  shift  of  maximum  cooling 
rate  away  from  surface  to  top  of 
cloud  bank  as  fog  bank  becomes 
radiatively  thick  (6  2. 

Oliver,  et  al.  (1978). 
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Figure  4.6.5  -  Turbulence  (q^) 
distribution  in  advective-radiative 
fog  of  Figures  4.6.3  and  4.6.4. 


4.7  Boundary  Layer  With  Stratus  Cloud 

A  boundary  layer  over  a  large  expanse  of  water  capped  by  a  subsidence  Inversion  will  turbulently 
Ingest  moisture  from  the  surface.  Long-wavelength  cooling  of  the  water-vapor-laden  air  may  then  lead  to 
condensation  somewhere  within  the  layer.  Calculations  using  relatively  similar  turbulence  closure  models 
have  been  made  by  Oliver,  et  al.  (1978),  Yamada  and  Mel  lor  (1980),  and  Moeng  and  Arakawa  (1980).  In  an 
expansive  basin  or  region  Influenced  only  by  advection  in  the  form  of  a  steady  subsidence,  a  periodic 
steady  state  may  form  in  which  a  solar  driven  diurnal  cycle  exists  about  a  mean  boundary  layer  state  with 
a  steady  mean  Inversion  height  controlled  by  the  radiative  conditions  (season  and  latitude),  surface 
temperature,  and  upper  conditions  (wind,  water  content  and  subsidence)  set  by  the  larger  circulation. 
Such  an  episode  was  Illustrated  by  Oliver,  et  al.  (1978)  for  conditions  approximating  the  simmer  stratus 
layer  off  the  coast  of  California. 
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The  conditions  for  this  illustration  were  selected  to  generally  conform  to  the  conditions  described 
by  Neiberger  (1944),  and  more  recently  observed  by  Mack,  et  al.  (1974).  Surface  boundary  conditions  are 
established  for  a  fixed  temperature  sea  state  of  17°C,  surface  water  mixing  ratio  at  saturation,  and 
surface  roughness  determined  by  Froude  scaling.  Radiation  conditions  are  for  latitude  40°N  at  simmer 
soltice.  Geostrophic  conditions  set  an  upper  level  wind,  Ug,  of  10  ms-1  and  subsidence  whose 
characteristic  value  is  0.5  x  10-5  s-1  (Neiberger,  et  al.  1961;  Lilly,  1968).  Initial  conditions  (which 
will  be  lost  after  several  days  of  simulation)  were  selected  to  correspond  to  a  mildly  stable  temperature 
profile  with  a  lapse  rate  suv/az  *  0.003°Cm-1  and  a  clear  sky  with  an  initial  relative  tumidity  of  0.9. 
The  calculation  is  begun  shortly  before  sunset  the  first  day  and  runs  for  six  days.  Behavior  during 
daylight  hours  must  be  considered  somewhat  approximate  in  this  illustration,  because  droplet  scattering  of 
direct  solar  radiation  which  would  reduce  solar  penetration  into  the  cloud  interior  was  neglected. 

The  quasi-periodic  evolution  of  cloudiness  in  the  boundary  layer  is  shown  in  Figure  4.7.1.  The 
stratus  grows  both  downward  and  upward  during  nocturnal  periods  and  thickens  until  sunrise,  when  the  solar 
heating  begins  to  disperse  the  stratus.  Radiative  cooling  is  dominant  at  the  cloud  top  and  produces  an 
unstable  lapse  rate  within  the  cloud.  A  correspondingly  enhanced  turbulence  production  is  also  maximum  in 
nocturnal  periods  and  in  the  upper  portion  of  the  cloud.  The  radiation  cooling  at  the  cloud  top 
strengthens  the  capping  inversion  which  oscillates  up  and  down  over  the  period  of  the  diurnal  cycle. 


suMwise 
Time,  hours 

Figure  4.7.1  -  Diurnal  variation  of  cloud  base  and  top  and  heat  and  moisture  fluxes  at 
surface.  Oliver,  et  al.  (1978). 


We  note  that  the  predicted  cloud  base  in  these  diurnal  cycles  is  highest  in  late  afternoon  and  lowest 
in  early  morning  (Figure  4.7.1),  while  cloud  top  is  highest  in  early  morning  and  lowest  in  late  afternoon. 
Correspondingly,  the  inversion  height  is  highest  in  early  morning  and  lowest  in  late  afternoon  -  a  result 
which  is  opposite  to  that  in  Section  4.1  for  inversion  height  cycles  driven  by  solar  heating  in  the 
cloud-free  boundary  layer.  This  predicted  result  is  in  accord  with  the  stratus  observations  of  Neiberger 
(1944)  as  well  as  those  of  Mack,  et  al.  (1974)  in  which  they  detected  a  diurnal  variation  of  inversion 
height  in  the  presence  of  California  coastal  stratus  which  regularly  showed  maxima  in  the  early  morning. 
It  is,  of  course,  the  radiative-turbulent  drive  of  cloud-top  cooling  which  allows  the  radiating  layer  to 
propagate  condensation  upward  during  the  night  while  turbulence  cools  the  cloud  Interior  below. 

Some  summary  aspects  of  this  illustration  in  comparison  with  the  observations  of  Mack,  et  al.  (1974) 
which  guided  the  choice  of  boundary  conditions  are  as  follows.  The  typical  liquid  water  content  of  the 
theoretically  generated  stratus  was  0. 1-0.2  g/kg  compared  to  0.1-0. 3  g/kg  observed.  The  average  height  of 
the  stratus  base  was  50-200  m  (theoretical)  compared  to  0-200  m  (observed).  The  average  height  of  stratus 
tops  at  maximum  inversion  were  600-800  m  both  predicted  and  observed.  The  theoretically  predicted  cooling 
rate  100  m  below  cloud  top  between  0100  and  0300  local  time  for  the  last  two  days  of  the  simulation  was 
0.0?3°Ch-1  compared  to  values  of  0.02b°Ch-1  between  0100  and  0500  observed.  The  amplitude  of  the  vertical 
diurnal  motion  of  the  cloud  base  was  100-150  m  predicted,  compared  to  observations  of  50-150  m. 

Me  have  seen  that  radiative-turbulent  coupling  of  subsidence  capped  stratus  layers  propagates  the 
base  downward  during  nocturnal  cooling.  If  the  base  propagates  completely  to  the  ground,  a  fog  is  formed 
at  the  surface.  This  class  of  fog  event  is  common  fn  coastal  California  during  periods  of  strong 
subsidence  and  has  been  described  by  Leipper  (1948)  and  Mack,  et  al.  (1974).  Several  atmospheric  fluid 
mechanisms  have  been  postulated  to  explain  this  class  of  fog  event  Including  vertical  convection  driven  by 
the  sea-breeze  land-breeze  cycle  (Neiberger,  1944),  and  radiative  cooling  of  stratus  during  the  night 
(Mack,  et  al.  1974).  The  single  necessary  condition  all  observers  note  for  this  class  of  fog  is  the 
existence  of  the  late  afternoon  or  early  evening  inversion  height  at  about  400  m  or  lower.  The  typical 
clear  summer  afternoon  temperature  profile  is  as  indicated  in  Figure  4.7.2,  as  measured  by  Mack, 
et  al.  (1974).  A  strong  upper  inversion  exists  above  the  inversion  height  zi  with  lapse  rate  yu  typically 
about  0.04'Cm-1.  A  transition  region  of  thickness  uz*40  m  exists  above  zi  with  lapse  rate  of  Q.O^Cm-1, 
while  below  the  inversion,  the  lapse  rate  is  -O.Q15<’Cm-1.  Figure  4.7.3  shows  the  resulting  liquid  water 
content  contours  of  a  calculation  which  begins  with  such  an  Initial  profile  over  a  fixed  sea  surface 
state,  as  calculated  by  Oliver,  et  al.  (1978).  Initial  relative  hunidity  is  0.97  in  the  region  up  to  the 
inversion  falling  off  to  a  value  corresponding  to  a  mixing  ratio  of  H  •  0.001  over  the  region  of  the 
transition  zone.  The  calculation  was  done  for  solar  inclination  at  solstice  and  begun  at  1900  local  time. 
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LIQUID  WATER  ISOPLETHS 


Figure  4.7.2  -  Typical  late 
temperature  profile  measured  by 
al.  (1974). 


afternoon 
Mack,  et 


Figure  4.7.3  -  Evolution  of  fog  resulting  from 
the  lowering  stratus.  Fog  shows  some  tendency 
to  form  at  surface  just  before  downward 
propagating  stratus  reaches  the  surface. 
01 iver ,  et  al .  (1978) . 


The  measurements  of  Mack,  et  al.  (1973)  of  liquid  water  content  in  stratus  lowered  to  the  surface  and 
for  conditions  typical  of  the  present  calculation  indicate  liquid  water  content  ranging  from  0.1  to 
0.2  g/kg  and  increasing  from  the  surface  upward  into  the  cloud.  This  trend  is  consistent  with  the 
calculated  results,  although  the  calculated  water  content  near  the  surface  is  less  than  measured  most 
likely  due  to  the  presence  of  drizzle  in  the  measured  case. 

There  is  a  large  difference  between  direct  solar  and  terrestrial  absorption  lenytbs.  Because  of  the 
strong  terrestrial  absorption  and  emission  of  the  water  droplets,  the  cloud  top  is  in  strong  radiative 
cooling  which  dominated  over  solar  heating  even  in  daylight  hours.  Direct  solar  radiation  is  absorbed 
much  deeper  within  the  cloud  interior.  Evaporation  by  solar  heating  occurs  rather  uniformly  throughout 
the  cloud  because  of  the  long  absorption  length  for  solar  radiation.  Strong  cloud-top  cooling  brought 
about  by  the  short  terrestrial  absorption  length  provides  the  instability  and  turbulence  production 
machinery  for  the  cloud  to  turbulently  transport  the  heat  deposited  by  the  sun  in  its  interior.  The 
details  of  this  turbulent  transfer  can  be  seen  in  Figure  4.7.4  where  the  turbulent  heat  flux  is  shown 
directed  from  cloud  base  to  cloud  top.  It  should  also  be  noted  that  heat  flows  turbulently  downward  above 
the  cloud  top  to  the  strong  radiatively  cooling  surface.  This  heat  flux  is  supported  by  the  turbulence 
which  has  diffused  a  few  tens  of  meters  above  the  cloud  top  into  the  strongly  stable  reyion  above  the 
cloud. 


»'•/,  m'C/bee 

Figure  4.7.4  -  Turbulent  virtual  heat  flux  for 
stratus  lowering  fog  at  time  of  stratus 
formation  (0130),  lowering  to  surface  (0500), 
and  lifted  state  (1500).  Oliver,  et  al.  (1978). 


It  is  noteworthy  that  the  corresponding  turbulent  moisture  flux  is  always  positive,  representing  a 
continuous  flow  of  water  from  the  sea  surface  upward  into  the  cloud.  At  no  time  (after  stratus  formation) 
during  the  course  of  lowering  and  raising  is  the  water  flux  negative,  indicating  a  diffusion  of  water 
downward  from  the  cloud  to  the  region  below  the  base.  Hence,  the  relative  downward  flow  of  water  from  a 
radiatively  cooling  cloud  could  only  take  place  by  precipitation  which  is  not  included  in  the  present 
calculation.  Precipitation  transport  of  water  down  to  the  subcloud  region  will  enhance  the  liquid  water 
level  in  that  region.  It  also  provides  a  mechanism  for  the  removal  of  liquid  water. 

The  calculations  provide  some  insight  as  to  the  appropriateness  of  the  entrainment  hypotheses  upon 
which  simple  cloud-top  mixing  models  are  based  (e.g.,  Lilly,  1968;  Schubert,  1976).  In  figure  4.7.4,  the 
distribution  of  the  turbulent  heat  flux  in  diurnally  varying  stratus  indicates  that  the  virtual  heat  flux 
has  a  positive  definite  average  over  the  layer  and  through  the  cloud  top  where  a  small  region  of  negative 
flux  exists.  Thus,  Ball's  hypothesis  that  the  average  should  vanish  (Ball,  1960)  and  the  corresponding 
maximum  entrainment  condition  which  follows  from  it,  are  not  descriptive  pf  a  cloudy  layer.  This 
turbulent-radiative  model  indicates  that,  on  average,  radiative  transport  is  as  important  as  turbulent 
transport  in  the  subcloud  region. 

The  overall  turbulence  levels  are  shown  in  Figure  4.7.6.  In  the  early  morning  when  stratus  covers 
only  a  small  portion  of  the  initial  mixed  layer,  the  turbulence  levels  have  a  typical  dry  mixed- layer 
structure.  After  the  cloud  base  rises,  the  confinement  of  radiative  cooling  to  the  cloud  top  and  the 
presence  of  heating  in  the  cloud  interior  skews  the  distribution  to  concentrate  and  locate  the  maximum 
turbulence  levels  near  the  cloud  top. 


Figure  4.7.5  -  Turbulence  levels  in  a  stratus 
lowering  fog  at  the  times  and  conditions  of 
Figure  4.  7.4. 


This  calculation  indicates  that  radiation  cooling  and  the  corresponding  turbulence  production  it 
effects  at  a  stratus  top,  cools  the  atmosphere,  and  propagates  the  cloud  top  upward  and  cloud  base 
downward.  The  effective  region  over  which  the  cloud  base  may  be  propagated  away  from  the  top  by 
turbulence  production  at  the  top  during  the  nocturnal  period  is  determined  by  the  long-wavelength 
radiation  absorption  length  and  by  the  turbulent  scale  length  over  which  turbulence  may  be  diffused 
downward  from  the  top  of  the  base.  For  these  typical  stratus  conditions  used  as  the  basis  of  the 
calculation,  we  find  the  region  to  be  850  m  (Figure  4.7.3).  The  empirical  observation  that  the  downward 
propagating  stratus  base  can  reach  the  ground  during  nocturnal  cooling,  while  the  cloud  top  is  at  400  m  or 
lower,  is  quite  consistent  with  this  result. 

Thermal  radiation  has  been  shown  to  play  a  significant  role  in  atmospheric  fog  and  stratur  cloud 
events.  Illustrations  of  the  dynamics  of  a  turbulent,  radiatively  coupled  low-level  atmosphere  :s-'  warm 
water  fog  and  subsidence-capped  stratus  show  consistent  quantitative  features  found  in  observatii-s  of 
such  fog  and  cloud.  Inversion  heights  with  cloud  in  the  boundary  layer  are  predicted  to  have  maxima  in 
early  morning  and  minima  in  late  afterr  i  In  contrast  to  the  diurnal  behavior  of  inversion  heights  in  a 
dry  boundary  layer  over  land.  Stratui  is  shown  to  propagate  both  downward  and  upward  during  nocturnal 
cooling  periods  under  the  Influence  of  the  radiant/turbulent  driving  mechanism  located  in  the  cloud  top. 
If  the  late  afternoon  inversion  height  is  low  enough,  the  stratus  Is  predicted  to  reach  the  surface. 
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V.  Organized  Features  within  the 
Atmospheric  Boundary  Layer 


In  even  quite  chaotic  turbulent  flow,  it  has  long  been  recognized  that  individual  organized  features 
may  be  discerned  if  the  flow  is  inspected  carefully  enough.  Leonardo  da  Vinci's  sketches  of  flow  in  the 
hydraulic  channels  of  his  day  clearly  show  this.  When  simulating  the  atmospheric  boundary  layer,  we  must 
decide  when  to  seek  to  resolve  particular  features  rather  than  to  ensemble  average  over  them.  Particular 
features  may  be  resolved,  using  a  model  of  the  type  discussed  in  Chapter  2,  by  incorporating  the  boundary 
conditions  for  the  feature  into  the  boundary  conditions  for  the  mean  variables. 

The  most  ambitious  calculations  of  this  type  have  been  Deardorff's  (1972,  1974a  and  b)  completely 
3-D,  unsteady  calculations  of  neutral  and  unstable  planetary  boundary  layers.  He  found  little  difference 
between  the  results  of  his  second-order  closure  calculation,  and  that  using  a  simpler  eddy  viscosity 
approach.  This  is  consistent  with  the  premise  that  he  was  able  to  use  a  sufficiently  fine  mesh  size  in 
his  calculation  so  that  the  major  contributions  to  the  turbulent  correlations  were  resolved  as  part  of  the 
mean  flow  calculation.  The  small  scale  motion  was  needed  to  dissipate  turbulent  energy,  but  the  details 
of  how  this  energy  was  dissipated  at  this  small  scale  appeared  relatively  unimportant  to  the  dynamics  of 
the  larger  scale  motion.  Our  purpose  herein  is  not  to  look  at  the  flow  in  this  great  detail,  but  rather 
to  consider  a  few  situations  where  organization  of  the  large  scale  flow  is  particularly  important. 

All  of  the  examples  considered  in  detail  in  this  chapter  involve  a  two-dimensional  symmetry,  either 
planar  or  ax i symmetric,  of  the  large  scale  flow.  The  three-dimensional,  smaller  scale  flow  must  be 
represented  by  the  closure  modeling.  Thus  we  expect  closure  to  be  less  critical  to  these  flows  than  it  is 
to  the  strictly  one-dimensional  simulations,  but  more  critical  than  it  is  to  a  fully  three-dimensional 
simulation.  The  particular  features  chosen  to  be  discussed  in  detail  reflect  a  personal  bias.  Other 
features  which  may  be  fruitfully  simulated  in  this  2-0  manner  are  the  coastal,  sea-breeze,  circulation 
(Pielke,  1974b;  Lewellen  and  Teske,  1976b);  frontal  dynamical  interaction  with  the  surface,  topography 
induced  features  (see  review  by  Hunt,  1980);  and  internal  wave  dynamics  in  stable  regions  of  the  boundary 
layer. 

Some  features  are  inherently  three  dimensional.  These  continue  to  present  a  real  challenge  to  the 
turbulent  transport  modeler.  Generally,  further  compromises  must  be  made  in  modeling  the  turbulence  in 
order  to  permit  the  added  nunerical  complexity  inherent  in  fully  3-D,  unsteady  calculations.  Yamada,  in 
his  3-D  simulations  of  flow  in  the  neighborhood  of  a  cooling  pond  (1979),  of  flow  in  a  complex  mountain 
terrain  (1980),  and  of  lake  breeze  circulations  around  Lake  Michigan  (1979),  has  used  the  one-dimensional 
turbulent  transport  model  given  in  Eqs.  (2.3.29  and  2.3.30).  Pielke  (1974a)  and  Tapp  and  White  (1976)  all 
use  O'Brien's  (1970)  eddy  viscosity  parameterization  in  their  3-D  simulations  of  the  sea  breezes  over 
south  Florida. 

The  brief  sampling  of  analyses  of  organized  boundary-layer  features  presented  within  this  chapter  are 
not  meant  to  imply  that  most  such  features  are  well  understood.  Quite  the  opposite  is  true.  Even  under 
conditions  typical  of  a  “horizontally  homogeneous"  fair  weather  boundary  layer  ,  Barnes,  et  al.  (1980) 
found  up  to  20%  of  the  "turbulent"  fluxes  to  be  transported  by  eddies  of  the  order  of  10  km  wide.  Much 
remains  to  be  learned  about  such  features. 


5.1.  Flow  Over  an  Abrupt  Change  in  Surface  Conditions 

The  simplest  type  of  feature  is  that  which  regularly  occurs  whenever  there  is  an  abrupt  change  in 
surface  conditions.  An  internal  boundary  is  induced  with  it's  thickness  growing  downwind  of  the  surface 
change.  Calculations  have  been  made  by  Rao,  et  al.  (1974)  and  by  Lewellen,  et  al.  (1974)  of  both  the 
transition  from  rough-to- smooth  and  smooth- to- rough  surfaces.  The  results  of  both  calculations  are 
compared,  in  Figure  5.1.1,  with  the  surface  data  of  Bradley  (1968),  who  recorded  velocities  and  surface 
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Figure  5.1.1  -  Response  of  surface  shear  to  step  change  in  surface  roughness.  Lewellen 

(1977)  (_ - );  Rao,  et  al.  (1974)  ( - ).  Data  from  Bradley  (1968).  (a)  Smooth  to  rough 

Zg  *  0.002  cm,  z0  *  0.25  cm;  (b)  rough  to  smooth  Zg  *  0.25  cm,  z0  *  0.002  cm. 
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shear  stress  at  a  little-used  airfield  in  Australia.  The  calculations  begin  with  the  velocity  and 
turbulence  components  in  equilibrium  for  one  value  of  aerodynamic  roughness,  Zq.  The  surface  boundary 
condition  is  changed  to  that  appropriate  for  the  new  value  of  z0  at  x  =  0.  When  both  runs  are  allowed  to 
continue  far  downstream,  the  turbulence  comes  into  equilibrium  with  the  new  value  of  20.  Rao,  et 
al.  (1974)  used  the  second-order  closure  model  developed  by  Wyngaard  and  Cote  (1974)  which  has  somewhat 
different  modeled  terms  than  those  given  in  Chapter  II  and  used  in  the  calculations  by  Lewellen,  et 
al.  (1974).  Both  calculations  show  reasonably  good  agreement  with  the  observations.  The  major  adjustment 
in  surface  shear  stress  is  made  within  a  few  meters,  but  a  distance  of  the  order  of  a  several  tens  of 
kilometers  is  required  for  the  complete  boundary  layer  to  make  its  transition. 

A  similar  sharp  change  in  the  boundaf y  layer  structure  occurs  when  there  is  an  abrupt  change  in 
surface  temperature.  This  is  often  seen  along  shorelines  where  temperature  differences  between  the  land 
and  water  may  easily  exceed  10°C.  This  leads  to  what  Lyons  (1975)  has  termed  a  thermal  internal  boundary 
layer  (TIBL).  Relatively  cold  air  flowing  over  a  warmer  surface  leads  to  a  growing  layer  of  increased 
turbulence  produced  by  buoyancy. 


5.2  Longitudinal  Roll  Vortices 

The  type  of  large  scale  coherent  structures  which  occur  most  often  in  the  boundary  layer  ar e  eddies 
which  extend  the  full  height  of  the  mixed  layer  in  the  form  of  longitudinal  roll  vortices  aligned  close  to 
the  direction  of  the  mean  wind.  An  example  of  the  clouds  resulting  from  this  motion  are  shown  in  Figure 
5.2.1.  When  the  boundary  layer  is  treated  as  spatially  homogeneous,  the  motion  of  these  roll  vortices 
must  appear  as  part  of  the  turbulence.  To  examine  the  relative  contributions  of  the  smaller  scale  random 
turbulence  and  the  relatively  organized  large  scale  rolls,  a  two-dimensional  calculation  was  performed  in 
the  (y.z)  plane  perpendicular  to  the  axis  of  the  rolls  (Lewellen,  et  al.  1980).  Periodic  boundary 
conditions  were  applied  so  that  the  domain  covers  a  wavelength  x  in  the  y  direction.  The  vertical  domain 
extended  from  the  surface  (represented  by  an  effective  roughness  zo)  to  slightly  above  the  capping 
inversion.  The  coordinate  system  permits  the  horizontal  roll  vortices  to  appear  as  part  of  the  ensemble 
mean  motion  when  calculating  the  unsteady  flow  in  the  unstable  boundary  layer. 


Figure  5.2.1  -  Apollo  photograph  of  Georgia  Coast, 
1700  Z  from  125  miles  high.  From  Brown  (1974). 
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The  nunerical  code  partitions  the  energy  between  the  mean  background  motion  which  is  a  function  of 
the  vertical  coordinate  only;  the  mean  quasi-periodic,  two-dimensional,  roll  vortex  motion  which  is  a 
function  of  z  and  y;  and  the  more  random  turbulent  motion  which,  although  three-dimensional  in  character, 
is  only  a  function  of  y  and  z  in  the  mean.  The  energy  in  the  organized  roll  motion  varies  with  the  ratio 
of  the  wavelength.  A,  of  the  roll  to  the  inversion  height,  z-j;  the  instability  of  the  layer  as  measured 
by  the  ratio  of  the  Monin-Obukhov  length,  L,  to  zj;  and  the  angle,  u,  between  the  roll  axis  and  the 
geostrophic  wind.  Figure  5.2.2  illustrates  the  cross-sectional  structure  of  the  stream  function  of  the 
roll  motion  for  an  angle  of  -10°,  A/zj  =  3  and  L/zj  «  -0.1  to  approximate  observations.  Surface 
conditions  are  chosen  appropriate  for  a  sea  surface  at  constant  temperature.  Results  are  scaled  in  terms 
of  the  characteristic  velocity 


to  minimize  the  influence  of  changes  in  wt)v.  (The  angled  bracket  represents  an  average  over  y.)  The 
geostrophic  velocity  is  held  at  *  10  m/s  throughout  the  calculations. 


indicate  flow  direction.  (Lewellen,  et  al.,  1980). 


The  corresponding  cross-sectional  structure  of  the  turbulent  kinetic  energy  and  the  tumidity  are 
shown  in  Figures  5.2.3  and  5.2.4.  The  mmerical  results  for  the  relative  contributions  of  the  mean  roll 
motion  and  the  roll  modulated  turbulence  to  the  transport  of  momentum  and  heat  compare  reasonably  (Teske 
and  lewellen,  1979)  with  the  observations  of  leMone  (1976). 

In  Figures  5.2.3  and  5.2.4,  the  background  levels  have  been  removed;  hence  aqz  shows  some  regions  of 
negative  values.  The  updraft  region  shpwn  in  Figure  5.2.2  produces  a  band  with  a  spread  of  about  0.25  A 
which  has  greater  than  background  C  and  qz.  The  rest  of  the  qz  region  is  dominated  by  a  much  broader  area 
of  less  than  average  values  of  aq*.  Figure  5.2.4  shows  a  noticeable,  almost  jet-like  character  to  aC. 
The  heat  flux  Is  concentrated  upward  in  the  same  region  also,  producing  a  temperature  overshoot. 

The  vertical  structure  of  the  flux  of  a  species  such  as  himidity  from  the  surface  is  shown  in  Figure 
5.2.5.  At  the  surface,  the  vertical  transport  is  all  carried  by  the  small  scale  turbulence  designated  by 
(wc),  but  at  z/zi  as  0.4,  approximately  2/3  of  the  vertical  flux  is  contained  in  the  mean  roll  motion 
designated  by  <WC>.  The  heat  flux  structure  is  similar  (lewellen,  et  al.  1980)  except  for  an  undershoot 
in  <We>  above  z/zi  as  0.75  forced  by  the  capping  inversion. 


Figure  S>.2.3  -  Perturbation  turbulence  contours  for  the  conditions  given  in  Figure  5.2.2 
Here  *  0.4  is  the  normalizing  maximum  value. 
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Figure  5.2  A  -  Roll  species  contours  for  the  conditions  of  Figure  5.2.2.  The  normalizing 
maximum  value  is  4C|Tvaxw./cw0  •  6.6<i. 


Figure  5.2.6  -  Composite  y-averaged  variation 
of  the  vertical  spegigs  flux  in  the  rolls  <WC> 
and  the  turbulence  <wc>  for  the  conditions  of 
Figure  6. 2. 2.  Here  cw0  is  the  species  flux  at 
V 


Figure  5.2.6  -  Comparison  betvieen 
the  one-dimensional  and 
two-dimensional  predictions  for 
vertical  velocity  variance.  For  the 
two-dimensional  calculation,  the 
y-averaged,  large  scale  roll 
contribution  is  denoted  by  <U^>  and 
the  smaller  scale  turbulence  by 
<ww>. 


If  the  closure  modeling  were  exact,  the  one-dimensional  (horizontally  homogeneous)  analysis  of  this 
problem,  which  averages  over  the  periodic  large  eddy  roll  structure  and  treats  the  roll  energy  as  part  of 
the  turbulent  kinetic  energy,  would  yield  the  same  answer  as  obtained  by  horizontally  averaging  over  the 
two-dimensional  result.  Since  the  tvo-dimensional  computation  allows  the  dominant  large  eddy  structure  to 
be  determined,  the  closure  modeling  should  be  less  critical  in  this  computation  than  it  is  in  the 
analogous  one-dimensional  computation.  In  fact,  a  2-D,  large  eddy  model  using  eddy  viscosity  by  Mason  and 
Sykes  (1980)  shows  many  of  the  same  features  demonstrated  here.  For  this  phenomenon,  which  has  a  strong 
two-dimensional  character,  our  model  represents  an  intermediate  step  between  depending  completely  on 
closure  modeling  and  depending  on  sub-9rid  closure  modeling.  Thus,  the  two-dimensional  result  can  be  used 
to  check  strengths  and  weaknesses  of  the  one-dimensional  model. 

The  separate  components  of  the  velocity  variance  are  shown  in  Figures  5.2.6  through  b.2.&.  There  is 
a  reasonable  comparison  between  the  vertical  variance  predicted  by  the  one-dimensional  model  and  the 
spatial  average  of  the  two-dimensional  results.  But  the  other  two  components  show  significant 
differences. 


figure  5.2. 1  -  Comparison  between 
the  one-dimensional  and 
two-dimensional  predictions  for 
transverse  velocity  variance.  For 
the  two-dimensional  calculation,  the 
y-averaged,  large  scale  roll 
contribution  Is  denoted  by  <V*>  and 
the  smaller  scale  turbulence  by 
<vv>. 


Figure  5.2.8  -  Comparison  between 
the  one-dimensional  and 
two-dimensional  predictions  .for 
longitudinal  velocity  variance.  For 
the  two-dimensional  calculation,  the 
y-averaged,  large  scale  roll 
contribution  Is  denoted  by  and 
the  smaller  scale  turbulence  by 
^uu>. 
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The  analogous  comparisons  of  the  vertical  momentum  flux  and  the  vertical  heat  flux  are  shown  in 
Figures  5.2.9  and  5.2.10.  In  the  middle  of  the  mixed  layer,  a  major  portion  of  the  turbulent  transport  is 
carried  by  the  large  roll  structure.  Both  models  predict  a  heat  flux  which  is  counter  gradient  over 
approximately  one-half  of  the  mixed  layer  depth.  The  two-dimensional  model  predicts  a  larger  undershoot 
in  the  heat  flux  near  the  inversion.  This  undershoot  in  the  heat  flux  is  forced  by  the  large  eddies 
overshooting  their  thermal  equilibrium  position.  The  two-dimensional  model  can  resolve  this  motion  but 
the  one-dimensional  model  must  rely  on  the  third-order  modeled  terms  to  correctly  incorporate  the  results 
of  these  large  eddies.  The  two-dimensional  model  is  also  able  to  account  for  more  influence  of  the  wave 
motion  in  the  stable  region  than  is  represented  in  the  one-dimensional  model.  This  large  difference  in 
the  two  models  is  also  apparent  in  Figure  5.2.11  which  shows  the  temperature  variance.  By  far,  the 
largest  contribution  to  temperature  variance  at  the  top  of  the  mixed  layer  is  the  wave  motions  induced  by 
the  large  rolls.  However,  even  the  small  scale  turbulent  variance  is  approximately  a  factor  of  two  larger 
than  that  predicted  by  the  one-dimensional  model  at  ih\  =  1. 


Figure  5.2.9  -  Comparison  between  the 
one-dimensional  and  tv«-dimensional  predictions  for 
the  vertical  longitudinal  momentum  flux.  The 
y-averaged  roll  contribution  is  denoted  by  <UW>  and 
the  turbulent  transport  by  <uw>. 


Figure  5.2.10  -  Comparison  between  the 
one-dimensional  and  two-dimensional  predictions 
for  the  vertical  heat  flux.  The  y-averaged  roll 
contribution  is  denot si  by  <Wt>>  and  the 
turbulent  transport  by  <w»>. 
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Figure  5.2.11  -  Comparison  between  the 

one-dimensional  and  two-dimensional  predictions  for 
the  temperature  variance.  The  y-averaged  roll 
contribution  is  denoted  by  <o^>  and  the  turbulent 
contribution  by  <«*>. 


Figures  5.2.6  through  5.2.11  show  that  although  the  turbulent  transport  compares  quite  favorably  for 
these  two  models,  there  is  considerably  more  kinetic  energy  near  the  surface  and  near  the  inversion  in  the 
more  correct  two-dimensional  model.  The  deficiencies  of  the  one-dimensional  model  appear  to  reflect  the 
difficulty  of  modeling  with  an  anisotropic  turbulent  macroscale.  It  Is  clear  from  the  two-dimensional 
simulation  as  it  Is  from  field  observations,  that  at  both  the  top  and  bottom  of  the  domain,  the 
characteristic  horizontal  scale  of  the  turbulence  is  much  larger  than  the  characteristic  vertical  scale. 
It  appears  that  the  next  step  in  improving  the  model  calls  for  Incorporating  some  structural  shape  in  the 
scale  representation,  as  discussed  in  Section  3.3. 

The  roll  structure  shown  in  Figure  5.2.2  should  be  expected  to  be  a  function  of  zi/L.  The  exact 
dependence  remains  to  be  unraveled,  but  l  would  expect  It  to  lead  to  some  modification  in  the  aspect  ratio 
of  the  rolls  and  their  alignment  with  respect  to  the  mean  wind.  A  better  definition  of  this  dependency 
can  possibly  permit  satellite  pictures  of  this  commonly  visible  boundary  layer  feature  to  yield  useful 
information  about  the  concurrent  boundary  - 1 ayer  structure.  Figure  5.2.12  shows  the  winter  time 
obervations  of  roll  clouds  over  the  Bering  Sea.  As  reported  by  Walter  (1980),  such  pictures  may  be 
analyzed  to  obtain  quantitative  Information  about  the  boundary  layer  thickness  and  its  rate  of  growth. 


Figure  5.2.12  -  (a)  NOAA-4  satellite  photograph  of 
the  Bering  Sea  showing  the  Ice  edge  and  roll  clouds 
(2127  GMT  19  March  1975).  Walter  (1980). 
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(o)  Map  of  the  Bering  Sea  area. 


5.3.  The  Thunderstorm  Downdraft  end  Associated  Gust  Front 

If  convection  is  much  stronger  than  advection  locally,  then  the  dominant  eddy  may  appear  as  a 
combination  of  ring  vortices  about  the  same  vertical  axis.  Natural  cumulus  clouds  and  plumes  from  large 
heat  sources  are  examples  of  such  vertical  plumes  driven  by  buoyancy.  We  consider  here  an  example  of  such 
an  axisymmetric  eddy  driven  by  negative  buoyancy.  This  downdraft  leads  to  horizontal  divergence  in  the 
boundary  layer.  The  positively  buoyant  plume  leads  to  horizontal  convergence,  which  will  amplify  any 
ambient  vertical  vorticity.  Thus  strong  upward  plumes  may  be  expected  to  strongly  interact  with  the 
vertical  vorticity  component.  The  extreme  example  of  this  type  of  flow,  the  tornado,  will  be  considered 
in  Section  5.4. 

The  gust  front  caused  by  the  cold  outflow  emanating  from  a  thunderstorm  has  been  modeled  by  Mitchell 
(1975)  using  a  constant  eddy  viscosity  model  and  by  Teske  and  Lewellen  (1977)  using  an  axisymmetric 
Implementation  of  the  second-order  turbulence  model  given  in  Chapter  Z.  The  simulated  flowfield  is 
Illustrated  in  Figure  5.3.1.  The  outflow  from  the  thunderstorm  is  Idealized  as  a  cold  jet  of  temperature 
e(r)  Impinging  normal  to  the  ground,  released  at  a  height  with  vertical  velocity  W(r).  The 
temperature  defect  below  the  ambient  temperature  is  caused  by  evaporation  of  falling  rain  by  relatively 
dry  air  at  some  altitude  z  >  z^x* 
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axlsymnetrlcal  thunderstorm  downdraft  interacting 
with  the  surface  to  form  a  gust  front. 


The  primary  variables  in  this  simulation  are  the  temperature  defect  (below  ambient)  of  the  jet,  its 
diameter,  and  the  height  at  which  it  is  released.  The  larger  scale  updraft  within  which  the  downdraft  is 
embedded  plays  a  role  in  retarding  the  late-time  gust  front  within  the  computational  domain,  but  should 
not  be  a  critical  factor  as  long  as  the  inflow  velocity  is  smaller  than  the  resulting  internal  simulated 
velocities.  The  other  critical  parameters  are  the  surface  temperature,  the  surface  roughness,  and  the 
stability  of  the  ambient  atmosphere. 

One  way  to  follow  the  structure  of  the  developing  gust  is  to  observe  the  movement  of  a  temperature 
isopleth,  in  this  case  Q  ■  -2°C,  as  shown  in  Figure  5.3.2,  for  the  simulation  sketched  in  Figure  5.3.1. 
At  t  =  0,  the  initial  linear  profile  is  shown,  but  by  t  *  200  sec,  the  structure  has  developed  into  a 
moving  front.  At  t  *  340  sec,  the  front  has  torn  away  from  the  downdraft  region.  At  later  times,  the 
strength  of  the  gust  decreases  as  the  effect  of  the  area  change  becomes  more  important.  However,  the 
height  of  the  gust  appears  to  grow  slowly,  so  that  by  t  a  1000  sec,  the  a  =  -2°C  line  reaches  nearly  1  Km 
in  altitude  at  a  distance  of  nearly  8  Km  from  the  source  centerline.  The  front  is  still  quite  strong  at 
this  point,  and  is  moving  outward  at  approximately  2.6  m/sec  with  maximum  horizontal  winds  within  the 
outflow  of  approximately  25  m/sec. 
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Figure  5.3.2  -  The  evaluation  of  a  simulated  gust  front  as  shown  by  isopleths  of  constant 
temperature  defect  (e  =  -2°C)  at  several  times  after  initialization.  (Teske  and  Lewellen, 
1977). 


When  this  numerically  simulated  front  is  visualized  in  terms  of  the  temperature  defect  intensity,  the 
leading  edge  appears  as  shown  in  Figure  5.3.3.  There  is  a  very  strong  qualitative  similarity  between  this 
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Figure  5.3.3  -  The  leading  edge  of  the  gust  front  shown 
in  terms  of  temperature  defect  intensity. 


in  Figure  5.3.2  when  visualized 


picture  and  the  gust  front  as  visualized  by  dust  in  Figure  5.3.4  (also  see  Turner,  1973).  The 
quantitative  predictions  of  maximum  mean  velocity,  velocity  shear  and  variance  also  appear  consistent  with 
available  observations.  As  part  of  a  study  to  evaluate  hazardous  conditions  for  aircraft  operations,  a 
sensitivity  analysis  has  been  performed  on  the  gust  front  as  a  function  of  five  different  dimensionless 
physical  parameters  (Teske  and  Lewellen,  1978).  We  will  return  to  this  in  Section  6.2. 


Figure  5.3.4  -  The  leading  edge  of  a  thunderstorm 
gust  front  as  marked  by  dust.  (Photo  by  Andrew 
Watson,  Miami,  Florida,  1975). 


5.4.  Tornado  Boundary  Layer 

The  horizontal  convergence  induced  by  positively  buoyant  plumes  or  updrafts  generally  leads  to  some 
amplification  of  ambient  vorticity.  At  the  barely  discernable  level  this  leads  to  the  small  whirls  in  the 
steam  fog  rising  from  a  warm  water  pond  on  a  cold  morning.  The  common  dust  devil  provides  a  somewhat 
larger,  and  more  impressive  demonstration  of  this  phenomenon.  Generally  larger  and  still  more  impressive 
is  the  water  spout.  Finally,  the  most  awesome  example  of  this  type  of  flow  is  the  mature  tornado. 

Recent  reviews  of  tornado  dynamics  have  been  given  by  Davies-Jones  (1980)  and  Lewellen  and  Sheng 
(1980b).  Although  the  precise  conditions  necessary  for  the  birth  of  a  tornado  are  still  subject  to  some 
speculation,  the  two  essential  ingredients  are:  (1)  an  ambient  conditional  instability  generated  by 
abundant  moisture  in  the  boundary  layer  capped  by  a  dry  air  mass  in  the  mid  troposphere;  and  (2)  ambient 
vertical  vorticity  in  the  form  of  a  mesocyclone  which  may  be  generated  by  the  amplification  of  the 
vorticity  associated  with  the  earth's  rotation.  Numerical  simulations  (e.g.,  Klemp  and  Wilhelmson,  1978) 
and  oual  Doppler  radar  observations  (e.g.,  Brandes,  1977)  also  demonstrate  a  critical  role  for  vertical 
shear  of  the  ambient  horizontal  winds  in  organizing  the  most  severe  storms. 

The  ambient  conditional  instability  is  released  through  condensational  heating  in  the  uplifting  of 
the  moist  surface  layer  air  and  evaporative  cooling  in  the  sinking  of  the  dry  mid-tropospheric  air.  The 
vertical  shear  in  the  horizontal  wind  allows  both  the  updraft  and  downdraft  to  become  well  developed 
before  they  interfere  with  each  other.  Both  the  Doppler  observations  and  the  numerical  simulations  also 
suggest  that  the  tilting  of  the  ambient  horizontal  vorticity  is  as  important  in  the  generation  of 
vorticity  as  the  convergence  term. 

A  vertical  jet  or  plume  issuing  into  a  horizontal  crossflow  generates  a  pair  of  vortices  with  equal 
magnitude  and  opposite  sign  circulation.  The  magnitude  of  the  vorticity  generated  in  this  manner  should 
be  of  the  order  of  the  product  of  the  distance  between  the  vortex  pair  and  the  relative  velocity  between 
the  plume  and  the  crossflow  divided  by  the  square  of  the  radius  of  the  region  of  concentrated  vorticity. 
Thus  a  10  m/sec  relative  velocity,  a  10  Km  distance  between  the  pair,  and  a  3  Km  region  of  concentrated 
vorticity  lead  to  a  value  of  vorticity  of  the  same  order  as  that  found  in  the  mesocyclone  (~  lO-^sec*1). 
As  long  as  the  vorticity  generated  by  this  tilting  process  is  of  the  same  order  as  that  resulting  from  the 
convergence  of  the  ambient  vorticity  due  to  the  earth's  rotation,  the  cyclonic  member  of  the  vortex  pair 
is  essentially  doubled  while  the  antlcyclonic  member  Is  nearly  eliminated.  The  remaining  single 
mesocyclone  tends  to  bear  to  the  right  of  the  ambient  wind  direction. 

A  fully  3-0  model  is  required  to  simulate  the  total  storm.  Including  all  the  Important  vorticity 
generating  terms.  However,  a  restricted  region  about  the  central  core  of  the  tornado  may  appropriately  be 
modeled  as  axisymmetric,  with  flow  into  and  out  of  the  boundaries  of  the  computational  domain.  An 
axisymmetric  version  of  the  basic  turbulent  transport  model  presented  in  Chapter  2  has  been  used  to 
simulate  this  tornado  core  boundary  layer  (Lewellen  and  Teske,  1977).  Figures  5.4.1  through  5.4.3  show 
the  results  of  a  computer  simulation  chosen  to  show  a  resemblance  to  the  1974  Xenia,  Ohio  tornado.  For 
this  simulation,  the  domain's  outer  radius  and  top  are  both  placed  at  400  m.  The  surface  roughness  is 
taken  to  be  large  (0.4  m)  to  correspond  to  a  relatively  urban  area.  All  components  of  the  velocity  are 
specified  at  the  outer  radius  Inflow  boundary  and  relatively  free  zero  slope  conditions  are  applied  at  the 
top  of  the  domain. 
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Figure  5.4.1  -  Mean  tangential  velocity  contours  as 
predicted  by  an  axi symmetric  tornado  model.  Contours 
labeled  in  tenths  of  maximum  V,  and  radial  and 
vertical  coordinates  normalized  by  the  radius  at 
which  Vmx  occurs. 


Figure  5.4.2  -  Mean  radial  velocity  contours  as 
predicted  by  an  axi symmetric  tornado  mode).  Contours 
labeled  in  tenths  of  maximum  V,  and  radial  and 
vertical  coordinates  normalized  by  the  radius  at 
which  occurs. 
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Figure  5.4.3  -  Mean  vertical  velocity  isopleths  for 
the  slnulated  tornado.  Contour  notation  is  the  same 
as  In  Figure  5.4.1,  except  normalization  is  by 
0.510  VMX. 
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The  mean  tangential,  radial,  and  vertical  velocity  fields  given  in  Figures  5.4.1  through  5.4.3  are 
normalized  by  the  maximum  tangential  velocity  which  occurs  at  r  s  180  m,  z  »  100  m.  The  absolute  value  of 
the  maximum  velocity  depends  strongly  on  the  relatively  uncertain  inflow  and  outflow  boundary  conditions. 
However,  when  the  model  result  is  scaled  to  give  a  Vwx  in  agreement  with  the  photogrammetric  analyses  of 
Umenhofer  and  Fujita  (1977),  the  vertical  distribution  of  tangential  velocity  at  r  *  180  m  is  given  in 
Figure  5.4.4.  Although  the  wind  scale  is  set  by  the  observations,  the  vertical  distribution  and  the 
fluctuations  about  the  mean  are  predicted  by  the  model.  Figure  5.4.5  shows  a  similar  representation  of 
the  vertical  velocity.  With  the  velocity  scale  set  in  the  previous  figure  there  is  no  arbitrariness  in 
the  predictions  of  Figure  5.4.5. 


Figure  5.4.4  -  Comparison  between  the  swirl  velocity 
profile  predicted  by  A.R.A.P.'s  model  at  r  =  180  m 
(Lewellen  and  Teske,  1977)  and  the  observations  of 
Umenhofer  and  Fujita  (1977).  Dashed  line  represents 
the  predicted  r.m.s.  fluctuation  about  the  mean  solid 
curve. 


Figure  5.4.5  -  Vertical  velocity  distribution  at  the 
same  radius  and  for  the  same  conditions  as 
Figure  5.4.4. 


Model  sensitivity  analyses  show  the  flowfleld  to  be  most  dependent  on  the  ambient  vertical  vorticity 
and  horizontal  convergence  occurring  in  the  parent  thunderstorm,  and  on  the  surface  roughness.  When  dual 
Ooppler  observations  (Brandes,  1979)  are  used  to  impose  the  boundary  values  of  vorticity  and  convergence 
at  1  Km  radius  and  1  Km  height,  the  model  Is  consistent  with  the  estimate  that  velocities  within  a  tornado 
will  very  rarely.  If  ever,  exceed  125  m/sec  (Lewellen  and  Sheng,  1980). 

A  dominant  feature  of  all  of  the  flowfields  computed  (lewellen,  et  a  1 . ,  1979;  and  Lewellen  and 
Sheng,  1979,  1980)  Is  that  the  maximum  winds  occur  at  quite  low  altitudes.  The  surface  layer  plays  the 
rather  paradoxical  role  of  increasing  the  velocity  in  the  neighborhood  of  the  ground  at  small  radii, 
before  actually  reducing  the  velocity  to  zero  at  the  surface.  The  strong  radial  inflow  in  the  surface 
layer  permits  the  streamlines  to  penetrate  to  smaller  radii  here  than  at  higher  altitudes.  Even  without 
total  conservation  of  angular  momentum,  this  permits  higher  swirling  velocities  to  develop.  The 
simulation  also  predicts  that  the  r.m.s.  average  of  the  velocity  fluctuations  reaches  values  as  high  as 
^*3  Vmax*  The  time  scale  of  the  eddies  is  such  that  the  average  strong  eddy  should  have  a  life  span 
somewhat  less  than  one  period  of  revolution  for  the  tornado.  This  large  magnitude  of  the  fluctuating 
velocities  indicates  that  they  should  be  considered  in  setting  design  criteria.  The  maximum  damage  will 
occur  in  swaths  where  the  instantaneous  fluctuating  velocities  add  to  the  local  mean  velocity.  Such 
damage  patterns  are  often  observed  (Fujita,  et  al.,  1976). 
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VI.  Example  Applications  of  Particular  Interest 


6.1.  Dispersal  of  Industrial  Plumes 

As  noted  in  the  introduction,  one  of  the  vital  services  that  the  atmosphere  provides  is  the  dispersal 
of  waste  products  from  the  machines  man  has  invented  to  serve  himself.  The  capacity  of  the  low-level 
atmosphere  to  perform  this  function  varies  strongly  with  meteorological  and  specific  site  conditions.  A 
large  number  of  semi -empirical  models  have  been  developed  to  predict  this  capacity  as  part  of  the  process 
of  deciding  what  facilities  may  be  built  where.  The  status  of  the  most  standard  of  these  models  have  been 
reviewed  by  Turner  (1980).  Our  purpose  herein  is  not  to  review  the  excellent  expositions  given  on 
atmospheric  diffusion  by  Pasquill  (1974),  Csanady  (1973),  and  Dobbins  (1979),  but  rather  to  briefly 
discuss  the  subject  in  the  context  of  turbulence  modeling. 

As  must  be  evident  from  the  preceeding  chapters,  the  motion  of  the  lower  atmosphere  has  a  rich 
variety  of  modes  for  transporting  and  intermixing  species.  Practical  methods  of  predicting  dispersal 
necessarily  aim  at  greatly  simplifying  these  mixing  processes  rather  than  attempt  to  solve  for  the  motion 
in  minute  detail.  In  attempting  to  describe  an  effluent  plume  of  a  nonreacting  conserved  species,  the 
most  important  features  are  the  mean  wind  which  determines  where  the  plume  will  go,  and  the  rate  of  spread 
of  the  plume  in  the  plane  perpendicular  to  its  mean  trajectory.  The  spread  of  the  plume  may  be  measured 
by 
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where  the  integrals  are  carried  out  over  the  plane  perpendicular  to  the  mean  trajectory  of  the  plume.  An 
equation  for  the  rate  of  spread  may  be  obtained  by  multiplying  Eq.  (8.1.9)  by  r2  and  integrating,  i.e.. 


Oof  1  D  <Cr2>  of  D  <C> 

Dt  "  <C>  Dt  ’  <C>  Dt 
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When  the  flux  goes  to  zero  on  all  the  boundaries  of  the  region,  this  last  equation  reduces  to 


From  Eq.  (2.3.21) 
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(6.1.5) 


and  when  <xju^c>  is  eliminated  between  Eqs.  (6.1.4  and  6.1.5)  a  second-order  differential  equation  for  o£ 
is  obtained 
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It  Is  instructive  to  look  at  Eq.  (6.1.6)  under  different  limiting  conditions.  First  consider  the 
simplest  limiting  condition  when  the  mean  wind  and  turbulence  quantities  are  homogeneous,  and  the 
gravitational  term  is  negligible.  In  this  limit 


D2of  ^  ZujUj  y  aC\  Aq  Dof 
Dt2  ’  '  <C>  X1  *Xj/‘  TUT 


(6.1.7) 
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which  is  satisfied  by  the  expression 


°r  “  Z^u^t  (6a-8) 


since 


<*i  -  -  «ij  <c>  (6.1.9) 


Equation  (6.1.8)  corresponds  to  the  classic  constant  diffusivity  solution  with  K  *  (4/3)Aq  and  the  plume 
spreading  as  t1/i!. 

The  initial  time  behavior  of  dispersion  from  a  point  source  does  not  generally  agree  with  that  given 
by  Eq.  (6.1.8)  even  when  the  ambient  atmospheric  turbulence  is  stationary  and  homogeneous. 
Deardorff  (1978)  noted  that  there  is  an  explicit  dependence  of  A/q  on  time  after  release  whether  closure 
is  made  at  first,  second,  or  third  order.  In  order  for  the  q  and  a  in  the  last  term  in  Eq.  (6.1.5)  to  be 
the  same  characteristic  velocity  and  length  scale  as  utilized  for  the  ambient  atmospheric  turbulence,  it 
is  necessary  for  the  plume  spread  to  be  larger  than  the  ambient  atmospheric  turbulent  length  scale.  This 
is  generally  not  the  case  in  the  initial  stages  of  a  release.  When  the  plume  spread  is  less  than  the 
length  scale  of  the  most  energetic  turbulent  eddies,  then  both  the  velocity  and  length  scale  of  the  plume 
may  be  significantly  smaller  than  the  comparison  scales  for  the  turbulent  velocity  field.  Lewellen  and 
Teske's  (1976)  approach  to  this  difficulty  was  to  take  the  plume  turbulent  scale  proportional  to  the 
spread  of  the  species  plume 


Ap  =  da  (6.1.10) 

subject  to  the  constraint 

Ap  <  At  (6.1.11) 


The  proportionality  factor  d  is  a  model  input  coefficient  of  the  order  of  1.  To  set  the  plume  velocity 
scale  qp,  they  assumed  that  most  of  the  interval  between  a^  and  Ap  would  fall  within  the  inertial  subrange 
of  the  turbulent  kinetic  energy  spectrum  so  that  the  rate  of  energy  dissipation  is  constant  as  the  energy 
cascades  from  the  larger  eddies  to  the  smaller.  From  Eq.  (2.8.5)  energy  dissipation  is  modeled  as 


c 


(6.1.12) 


so  that  it  is  appropriate  to  take 


dp 


(6.1.13) 


Equation  (6.1.7)  can  then  be  written  as 
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(6.1.14) 


The  late-time  plume  is  appropriately  described  by  Eq.  (6.1.8),  but  in  the  initial  phase  when  Ap/A^  «  1, 
the  last  term  in  Eq.  (6.1.14)  can  be  neglected  and 


of  *  U.jU^  t* 


(6.1.15) 


This  corresponds  to  the  early  time  behavior  predicted  by  statistical  theory  (Pasquill,  1974).  Thus 
Eq.  (6.1.7)  provides  a  smooth  transition  from  hyperbolic,  wave-like  dispersion  of  the  plume  when  the  plume 
spread  is  very  small  compared  to  the  turbulent  scale,  to  simple  parabolic  diffusion  when  the  plume  is 
large  enough  to  contain  the  motion  of  the  most  energetic  eddies. 

The  influence  of  buoyancy  enters  dispersion  in  two  ways.  First  it  has  a  strong  effect  on  the 
turbulent  velocity  components.  But  as  seen  in  Eq.  (6.1.6),  it  also  has  a  separate  effect  depending  upon 
the  covariance  between  the  temperature  and  the  species  concentration. 
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To  exhibit  the  influence  of  the  buoyancy  term  more  clearly,  let  us  consider  dispersion  in  the  tree 
convection  layer  in  the  absence  of  any  mean  velocity.  Then  the  equation  for  the  vertical  flux  of  any 
species  Eq.  (2.3.4)  may  be  written  as 


In  the  absence  of  significant  diffusion  contributions  from  the  third  term  on  the  r.h.s.  of  Eq.  (6.1.16),  a 
zero  occurring  in  wc  in  coincidence  with  a  zero  in  3C/3 z  cannot  be  maintained  unless  the  buoyancy  term  is 
also  negligible.  Conversely,  it  is  possible  for  a  vertical  gradient  in  mean  concentration  to  be 
maintained  in  response  to  the  buoyancy  term  even  in  the  absence  of  any  vertical  species  flux.  If  the 
covariance  of  species  and  temperature  is  positive  then  Eq.  (6.1.16)  shows  that  a  positive  3C/3z  can  be 
maintained  while  a  negative  covariance  would  lead  to  a  negative  ac/az.  A  pulsed  surface  release  into  the 
free  convection  layer  should  yield  a  positive  cs  and  consequently  a  positive  aC/az  in  the  immediate 
vicinity  of  a  nonabsorbing  surface.  This  leads  to  the  position  of  the  maximum  mean  concentration  rising 
from  the  surface,  as  observed  in  Figure  6.1.1  experimentally  by  Deardorff  and  Willis  (1974),  and  in 
Figure  6.1.2  as  numerically  calculated  by  Lewellen  and  Teske  (1976).  Conversely  a  pulsed  release  near 
the  top  of  the  convectively  mixed  layer  will  yield  at  least  a  temporary  negative  co  and  the  position  of 
the  maximum  concentration  will  tend  to  fall.  Such  a  behavior  could  not  be  predicted  by  a  two-equation 
turbulence  model  as  long  as 


wc 
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Figure  6.1.1  -  Centerline  concentration  isopleths 
downwind  of  a  continuous  point-source  release  into  a 
free  convective  layer;  the  results  of  Experiment  A 
of  Deardorff  and  Willis  (1974). 


Figure  6.1.2  -  Ensemble  average  of  the  centerline 
concentration  Isopleths  downwind  of  a  continuous 
point-source  release  into  a  free  convective  layer,  as 
calculated  by  Lewellen  and  Teske  (1976). 


The  popular  Gaussian  plume  models  described  by  Turner  (I960)  depend  upon  an  empirical  solution  to 
Eq.  (6.1.6).  The  spread  o  is  generally  tabulated  as  a  varying  power-law  of  t  (or  x)  for  five  different 
stability  classes.  Some  of  the  better  methods  provide  for  allowing  at  least  the  horizontal  variance  to  be 
directly  measured  at  the  site.  Equation  (6.1.6)  appears  to  potentially  provide  the  basis  for  using 
turbulence  data  measured  at  the  specific  site  of  a  release  for  predicting  the  dispersal  of  that  plume 
downwind  of  the  site. 

Sample  second-order  closure  calculations  of  dispersion  within  the  atmospheric  boundary  layer  have 
been  made  by  Lewellen  and  Teske  (1975,  1976a  and  b),  Tamada  (1979a),  and  Teske,  et  al.  (1978),  but  no 
rigorous  comparisons  with  data  have  been  made.  Comparisons  between  calculations  and  observations  of  real 
field  data  are  compounded  by  the  fact  that  the  model  provides  a  simulation  of  the  ensemble  average  plume 
while  the  field  observation  provides  one  specific  realization  from  a  number  of  different  plumes  which 
might  be  members  of  the  simulated  ensemble.  These  two  should  only  be  expected  to  coincide  if  (1)  the 
meteorology  is  sufficiently  steady  for  the  turbulence  to  be  stationary,  and  (2)  the  time  over  which  the 
plume  is  averaged  is  long  compared  with  the  characteristic  time  of  the  large  turbulent  eddies.  Such  a 
long  time  average  over  a  stationary  plume  should  approximate  the  ensemble  average  of  a  number  of 
individual  plumes,  each  of  which  satisfy  the  same  initial  mean  values.  This  combination  provides  a 
stringent  test  for  any  data  set.  As  seen  in  the  earlier  chapters,  the  first  condition  severely  limits  the 
upper  bound  on  the  observation  time,  while  the  second  pushes  up  the  lower  bound  on  the  observation  time. 

As  an  example,  consider  the  turbulence  time  scales  likely  to  occur  under  unstable  afternoon 
conditions.  If  the  mixed  layer  height  is  of  the  order  of  1  km,  and  w.  is  of  the  order  of  1  m/sec,  then 
the  time  period  for  the  rotation  of  such  an  eddy  will  be  of  the  order  of  1  hr.  A  time  average  of  three 
hours  or  more  would  be  required  to  provide  the  equivalent  of  an  ensemble  average.  Samples  averaged  over 
times  of  one  hour  or  less  should  be  expected  to  show  considerable  variation  from  the  ensemble  mean  value. 
Some  indication  of  this  difference  can  be  obtained  by  calculating  the  plume  concentration  variances  from 
Eq.  (2.3.22) 


ac^  ac*  _  aC  a  /  ac^\  2bsqc* 


(2.3.22) 


Numerous  efforts,  currently  underway,  to  evaluate  the  accuracy  of  existing  dispersion  models 
(e.g..  Hi  1st,  1978;  Hi  Iyer,  et  al.  1979;  Ruff,  et  al.  1979;  Londergan,  1980;  Fox,  1981;  and  Pepper, 
1981)  may  lead  to  the  incorporation  of  some  type  of  variance  estimates  in  a  new  generation  of  dispersion 
models. 


6.2.  Hind  Shear  Around  A’rports 

Wind  shear  ana  turbulence,  oa.  ticul  arly  at  the  lower  altitudes,  has  been  identified  as,  at  least,  a 
contributing  factor  tc  a  > jmber  of  different  aircraft  accidents  occuring  during  either  takeoff  or  landing. 
Flowfield  models  such  as  these  uiscussed  in  the  preceeding  five  chapters  are  useful  in  evaluating  the 
potential  hazard  preaentad  oy  wind  shear  and  turbulence  under  different  conditions.  A  second-order 
closure  model  of  the  atmospheric  boundary  layer  has  been  used  to  define  the  bounds  on  low-level  wind  shear 
under  different  conditions  (Lewellen  and  Williamson,  1976).  and  to  investigate  the  conditions  existing  at 
the  time  of  a  number  of  accidents  (Williamson,  et  al.  (1977). 

Even  within  a  one-dimensional  approximation  there  are  a  rich  assortment  of  possible  wind  and 
turbulence  profiles.  A  rough  summary  of  the  turbulence  intensity  for  different  conditions  is  shown  in 
Figure  6.2.1.  Refer  to  Table  6.2.1  for  the  key  to  the  different  conditions.  The  turbulence  is  normalized 
by  the  maximum  wind  velocity  between  the  ground  and  1  km  and  plotted  as  a  function  of  stability.  As 
indicated  in  Chapter  3,  it  is  not  possible  to  reduce  the  detailed  influences  of  stability  to  a  single 
parameter.  But  for  purposes  of  indicating  only  the  rough  dependence  of  stability,  we  have  chosen  a  single 
parameter  which  does  partially  combine  the  influence  of  the  surface  heat  flux,  which  controls  stability  in 
the  surface  layer,  and  the  influence  of  the  upper-level  potential  temperature  gradient,  which  influences 
the  elevation  of  the  inversion  layer  which  caps  the  boundary  layer.  This  parameter  is  a  bulk  Richardson 
number  based  on  the  velocity  and  potential  temperature  differences  between  the  surface  and  1  km.  The 
hysteresis  loops  presented  for  the  different  periodic  cycles  considered  indicate  the  type  of  scatter  to  be 
expected  in  attempting  to  predict  the  fluctuations  in  this  manner.  As  the  bulk  Richardson  number 
approaches  -«,  corresponding  to  free  convection,  the  normalized  q  will  approach  -.  On  the  stable  side, 
increasing  Ri ( 1000)  tends  to  decrease  the  normalized  q  whether  the  increase  in  stability  is  caused  by 
increased  surface  cooling  or  a  decrease  in  the  upper-level  inversion  height. 

A  similar  rough  summary  of  the  maximum  wind  shear  occurring  at  or  above  10  m  altitude  in  the  boundary 
layer  is  given  in  Figure  6.2.2.  There  is  a  fairly  consistent  trend  of  increasing  shear  with  increasing 
stability,  although  there  is  a  relatively  strong  Influence  of  thermal  wind  which  tends  to  mask  this  trend 
to  some  extent.  The  normal ization  of  q  by  the  maximum  velocity  within  the  layer  tends  to  roughly 
compensate  for  the  thermal  wind  influence  in  Figure  6.2.1.  Of  course,  the  wind  shear  in  a  profile  at  a 
given  time  cannot  be  characterized  by  a  single  number.  The  actual  changes  in  airspeed  which  an  aircraft 
will  encounter,  will  also  depend  on  the  vertical  depth  over  which  this  large  shear  exists.  This  must  be 
determined  by  looking  at  the  detailed  profiles.  Figure  6.2.2  does  indicate  the  rough  rule  of  thumb  that 
as  either  horizontal  temperature  gradients  Increase  or  vertical,  stable,  temperature  gradients  Increase, 
problems  with  wind  shear  may  be  expected  to  Increase. 
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Table  6.2.1.  Key  for  summary  plots  (Figures  6.2.1  and  6.2.2). 


1  For  the  conditions  corresponding  to  Figures  4.1.2  and  4.1.3 

4  For  the  conditions  corresponding  to  Figures  4.2.1  and  4.2.2 

5  Same  as  4  but  Mq  increased  from  1.25  to  10 

6  Same  as  4  but  upper  level  lapse  rate  increased  from  3°C/km 


to  6°C/km 


7 

Neutral  steady  state 

(R0  =  10s) 

8 

Neutral  steady  state 

<R0  =  10b) 

9 

Neutral  steady  state 

(R0  -  10?) 

10 

Stable  atmosphere 

(R0  *  10s. 

Ri  (10)  *  0.014) 

11 

Stable  atmosphere 

(R0  =  10s, 

Ri(10)  *  0.0036) 

12 

Stable  atmosphere 

(R0  ■  10s. 

Ri (10)  s  0.047) 

13 

Unstable  atmosphere 

(R0  -  10s. 

Ri(10)  =  -0.013,  Z|f/U, 

14 

Unstable  atmosphere 

r> 

o 

r4 

■  1 

o 

cz 

Ri  (10)  *  -0.U35,  Zi/Ug 

■  0.22) 

0.17) 


In  general,  the  wind  shear  decreases  with  altitude  in  the  lower  portion  of  the  boundary  layer.  This 
is  shown  in  Figure  6.2.3  to  supplement  the  summary  information  in  Figure  6.2.2.  The  solid  curve  shows  the 
distribution  of  wind  shear  with  altitude  for  the  neutral,  barotropic  distribution  of  Figure  3.2.1.  The 
bars  represent  the  change  in  wind  shear  caused  by  two  orders  of  magnitude  change  in  Ro.  The  dashed  curves 
represent  the  bounds  in  the  variation  of  wind  shear  caused  by  expected  values  of  Ri.  In  general, 
increasing  Ri  will  tend  to  increase  wind  shear,  as  seen  in  Figure  6.2.2,  but  the  minimum  wind  shear  need 
not  coincide  with  the  minimum  Ri  due  to  hysteresis  effects.  Barocl inicity  also  will,  in  general,  increase 
wind  shear.  It  is  not  possible  to  correlate  the  variations  indicated  in  Fiyure  6.2.3  with  time  of  day. 
Nocturnal  hours  do  tend  to  be  more  stable  (Ri  >  0)  and  afternoon  hours  more  unstable  (Ri  <  0),  but  the 
barocl inicity  is  associated  with  changing  weather  and  does  not  correlate  with  time  of  day.  A  companion 
curve  like  Figure  6.2.3  for  the  turbulence  intensity  as  a  function  of  altitude  would  be  too  misleading, 
since  there  is  no  consistent  trend  of  turbulence  intensity  with  altitude  over  parameter  space.  Under 
stable  conditions,  the  turbulence  will  decrease  with  altitude,  but  for  unstable  conditions  at  low  levels, 
it  will  increase  with  altitude. 


Time,  sec 

Figure  6.2.3  -  Summary  plot  of  wind  shear  as  a  function  of  altitude  showing  the  influence 
of  the  expected  range  of  variations  in  Ro,  Ri  and  rt. 


Many  cases  in  nature,  where  wind  shear  and  turbulence  are  a  problem  to  aviation,  are  further 
complicated  by  the  added  dimensionality  of  the  phenomena.  For  example,  the  strong  downdraft  associated 
with  flow  from  a  thunderstorm  discussed  in  Section  5.3  presents  a  horizontal  wind  shear  which  can  be 
extremely  hazardous.  A  slender  downdraft  can  persist  to  a  low  altitude  and  create  an  outflow  which  has 
horizontal  wind  speed  changes  in  excess  of  30  m/sec  within  a  few  kilometers.  Figure  6.2.4  shows  the 
relation  of  a  thunderstorm  downdraft,  as  deduced  by  Fujita  and  Caracena  (1977),  to  an  Eastern  Airlines 


65 


aircraft  on  Its  fatal  final  approach  to  JFK  June  24,  1975.  Analysis  of  the  flight  recorder  data  indicated 
that  the  combination  of  the  downdraft  and  the  sharp  change  In  horizontal  wind  speed  resulting  from  the 
outflow  from  the  downdraft  were  sufficient  to  cause  the  crash  approximately  1  km  short  of  the  runway.  The 
relative  contributions  of  the  negative  vertical  wlndspeed  and  the  sharp  dropoff  In  the  speed  of  the 
horizontal  head  wind  are  subject  to  question.  McCarthy,  et  al.  (1979),  argue  that  the  most  critical  wind 
shear  is  the  change  in  the  horizontal  wind  speed  on  a  time  scale  of  the  order  of  the  period  of  the  phugold 
oscillation  for  the  particular  aircraft. 
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Figure  6.2.4  -  Paths  of  14  aircraft  in  25  min  at  JFK  Airport  on  24  June,  1975,  as  analyzed  by  Fujita  and 
Caracena  (1977).  Each  path  was  shifted  toward  the  MNW  at  15  m  s*1  (30  kt)  to  convert  the  time  into  the 
space  relative  to  the  spearhead  echo.  The  echo,  as  seen  by  the  Atlantic  City  radar,  was  32  km  (20  mi) 
long  and  13  km  (8  mi)  wide,  covering  the  entire  area  between  LOM  (Localizer  Outer  Marker)  and  the  north 
end  of  runway  22-L.  There  were  downburst  cells  (DBCs)  along  the  south  edge  of  the  spearhead  echo.  Five 
aircraft  took  off  inside  the  sea  breeze  without  being  affected  by  downbursts. 

The  strongest  wind  shear  conditions  associated  with  a  thunderstorm  downdraft  will  occur  under  the 
center  of  a  relatively  small  diameter  downdraft.  Fujita  (1980)  calls  such  an  intense  narrow  downdraft  a 
"downburst'1  or  “mlcroburst".  The  downdraft  needs  to  be  sufficiently  large  that  the  core  can  penetrate  to 
the  surface  before  the  potential  temperature  defect  is  eroded  by  turbulent  mixing,  yet  sufficiently  narrow 
that  the  maximum  radial  outflow  occurs  at  a  small  radius.  Under  these  conditions  the  horizontal  wind 
shear  will  scale  as  [gA8/I0h] 1/2  where  At)  is  the  potential  temperature  defect  between  the  core  of  the 
downdraft  and  the  su  'ounding  ambient  air,  and  h  is  the  height  from  which  the  downdraft  originates. 

The  gust  front  which  may  emanate  for  some  distance  from  the  downdraft  presents  a  related,  but 
somewhat  different  wind  shear  hazard  (Goff,  1976).  A  detection  system  based  on  identifying  the  pressure 
signal  occurring  at  the  surface  under  such  a  discrete  feature  has  been  developed  by  Bedard,  et  al.  (1979). 
Figure  6.2.5  compares  a  pressure  trace  they  observed  at  Dulles  Airport  with  the  pressure  pulse  associated 


15  min 


Figure  6.2.5  -  Pressure  field  beneath  a  thunderstorm  downdraft: 

(a)  as  measured  on  26  June  1978  by  Bedard,  et  al.  (1979);  and 

(b)  as  given  in  numerical  simulation  by  Teske  and  Lewellen  (1978). 
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with  the  model  simulation  discussed  In  Section  S.3.  The  traces  are  quite  similar.  In  both  cases  there  is 
a  small  rise  followed  by  a  small  drop  preceedlng  the  main  pressure  pulse. 

The  results  of  model  studies  such  as  those  discussed  herein  can  provide  a  physically  consistent  set 
of  wind  and  turbulence  profiles  which  may  be  used  to  help  understand  what  meteorological  conditions  lead 
to  hazardous  low-level  wind  shear  and  turbulence  profiles,  and  to  provide  model  inputs  Into  flight 
simulators  for  pilot  training.  The  physical  constraints  provided  by  the  conservation  equations  should 
provide  a  valuable  supplement  to  simulator  models  based  on  correlations  of  empirical  data  (e.g..  Frost, 
Camp  and  Wang,  1978). 


6.3.  Refractive  Index  Fluctuations 

Turbulent  fluctuations  in  temperature  and  hunidity  force  fluctuations  in  the  atmospheric  refractive 
index  which  can  interfere  with  coherent  propagation  of  optical,  microwave,  or  acoustic  radiation.  These 
effects,  although  generally  detrimental  to  communication  systems,  have  formed  the  basis  for  remote  sensing 
of  the  atmosphere.  A  brief  review  of  the  rapidly  expanding  use  of  these  remote  sensors  in  boundary  layer 
meteorology  is  given  by  McBean,  et  al.  (1979).  Engineers  dealing  with  either  their  beneficial  results  or 
their  interference  effects  need  methods  for  calculating  these  fluctuations  in  the  atmospheric  boundary 
1 ayer. 

The  reflectivity  of  a  distributed  source  is  proportional  to  the  refractive  index  structure  function 
defined  as  (Tatarskii,  1971) 


C*  *  Cn(x)  -  n(x+  r)]2 

ri/t 


(6.3.1) 


where  n  is  the  index  of  refraction,  r  is  the  separation  distance  and  the  overbar  denotes  an  ensemble 
average.  For  many  applications  of  interest  in  the  atmospheric  boundary  layer  the  scattering  length  of  the 
propagating  wave  of  interest  will  fall  within  the  inertial  subrange  scale  of  the  turbulent  spectrun  Within 
this  subrange  where  essentially  isotropic  turbulence  exists,  C|(  will  be  independent  of  r. 

The  basic  turbulent  correlations  predicted  by  a  model  such  as  the  one  given  in  Section  2.3  can  be 
used  to  predict  the  refractive  index  structure  function,  C„.  The  combination  of  primary  variables, 
temperature  and  humidity  variances  and  covariances  required  to  compute  C*  depends  upon  the  wavelength 
being  refracted.  A  formula  for  cm-wavelength  radar  was  derived  by  lewellen  and  Teske  (1975) 
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where  B  »  7,730,  p  is  in  atmospheres,  T  is  temperature  in  °K  and  H  is  water  vapor  mixing  ratio.  The 
dominant  term  in  this  equation  for  moderately  moist  boundary  layers  is  proportional  to  the  hunldity 
fluctuations.  Fluctuations  in  the  mixing  ratio  of  the  order  of  a  fraction  of  a  gram/Kg  across  a  turbulent 
layer  of  the  order  of  ten  meters  thick  are  detectable  by  sensitive  radars. 

The  corresponding  fornula  for  optical  radiation  is 
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Expressions  for  acoustic  radiation,  as  well  as  optical  and  microwave,  have  been  given  by  Wesely  (19761. 
His  expressions  relate  cf  to  the  individual  temperature  and  water  vapor  pressure  structure  parameters,  cf, 
Cg,  and  Cej.  His  expression  for  radar  is  slightly  different  than  that  given  in  Eq.  (6.3.2)  since  he  does 
not  Include  the  influence  of  a  mean  humidity  level,  which  may  be  Important  in, the  marine  boundary  layer. 
Figure  6.3.1,  taken  from  Burk  (1981)  gives  the  results  of  a  calculation  of  cd  throughout  the  planetary 
boundary  layer  using  a  model  (Burk,  1977,  1980)  similar  to  the  one-dimensional  model  presented  in 
Chapters  II  through  IV. 


Figure  6.3.1  also  Includes  the  mixed  layer  scaling  variation  for  cd  given  by  Wyngaard  and 
Lemone  (1980)  along  with  their  data  collection.  They  also  provide  scaling  relations  for  the  structure 
functions  In  the  Immediate  neighborhood  of  the  inversion  based  on  Deardorff's  (1979)  idealization  of  the 
structure  of  the  interfacial  layer  defined  as  the  region  between  the  height  at  which  the  virtual  potential 
flux  first  goes  to  zero  and  the  height  at  which  all  the  turbulence  quantities  vanish.  These  expressions 
may  be  written  as: 
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„  ABv  Bv«  (6.3.6) 

bv  ~  Ziu  3 


2.2(aH)  by* 

Ctty,H  *  zZ/3 


where  aH  and  aov  are  the  changes  In  moisture  and  virtual  potential  temperature  across  the  inversion. 
These  relations  appear  to  depend  upon  the  dynamics  of  the  inversion  layer  being  strongly  coupled  to  the 
surface  dynamics  through  the  parameters  Zj  and  ev*.  This  need  not  always  be  the  case.  We  may  expect  tne 
turbulence  in  the  inversion  to  often  be  driven  by  local  internal  wave  breaking. 


Figure  6.3.1  -  Humidity  structure  function 
comparison  by  data  and  model  simulation  as  given  by 
Burk  (I960).  Data  collection  and  dashed  line,  the 
mixed  layer  scaling  prediction,  from  Wyngaard  and 
lemone  (1980).  Upper  scale  applies  only  to  solid 
line. 


Figure  6.3.2  gives  the  contours  of  CjJ  predicted  to  be  associated  with  the  simulation  given  in 
Section  5.2  for  one  set  of  conditions.  (Sea  surface  temperature  *  2U°C,  Ug  *  10  m/sec,  a  *  -10, 
i\! L  *  -10,  ae/iz  (above  mixed  layer)  *  0.01°/m,  H„  -  2.5  g/Kg.)  This  figure  shows  that  the  highest 
values  of  Cj|  occur  near  the  surface  and  In  the  vicinity  of  the  inversion  in  agreement  with  the 
one-dimensional  simulation,  but  there  Is  also  considerable  variation  across  the  roll.  The  large  values 
near  the  inversion  are  localized  in  the  region  of  the  rolls  where  the  updraft  of  moist  air  Interacts  most 
strongly  with  the  dry  air  above  the  Inversion.  Under  such  conditions  the  radar  return  can  show  the 
convective  streets  exemplified  In  Figure  6.3.3. 

The  combination  of  remote  sensors  capable  of  looking  in  great  detail  at  the  structure  within  the 
atmospheric  boundary  layer  and  models  able  to  calculate  the  turbulent  dynamics  seen  by  these  sensors 
should  greatly  Improve  our  understanding  of  the  rich  variety  of  fluid  dynamical  structures  which  occur 
particularly  in  the  neighborhood  of  the  top  of  the  boundary  layer. 
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Figure  6.3.2  -  Contours  of  the  radar  refractivity  structure  function  (C* ) 
svithin  the  roll  vortices  predicted  by  a  model  run  with 
T(sea  surface)  =  20°C,  Ug  =  10  m/sec,  Zj/L.«-10,  o  =  -10°, 
30/dz(above  inversion)  =  0.01  C/m,  Ah  (inversion)  =  2.5  g/Kg. 

Contours  labeled  n  denote  2*n~7l  x  10rlJ  so  that  there  is  a  factor 
of  2  difference  between  each  of  the  contours. 


Figure  6.3.3  -  PPI  photograph  taken  with  an  S-band 
radar  at  Wallops  Island  showing  alignment  of 
convective  cells  at  the  top  of  the  mixed  layer. 
(Konrad,  1970). 


6.4.  Ory  Deposition 

The  final  application  considered  here  Is  that  of  the  deposition  of  a  trace  species  from  the 
atmosphere  to  the  surface.  This  is  important  to  the  depletion  of  plumes  transported  over  long  distances. 
Deposition  is  the  negative  of  surface  evaporation,  a  phenomenon  of  utmost  Importance  to  agricultural 
meteorologists.  For  common  Industrial  pollutants,  approximately  half  of  the  pollutant  is  eventually 
transported  to  the  surface  directly  by  dry  deposition  with  the  rest  of  it  washed  out  of  the  air  by 
precipitation.  Thus  dry  deposition  Is  an  important  component  of  what  has  become  popularly  known  as  the 
acid  rain  problem. 
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A  number  of  reviews  of  deposition  have  been  made  (e.g.,  Sehmel,  1980;  SI  inn,  1976;  Chamberlain, 
1979;  and  McMahon  and  Denison,  1979).  In  spite  of  the  relatively  large  amount  of  data  accumulated,  the 
data  can  rarely  be  reliably  extrapolated  to  general  field  conditions,  due  either  to  insufficient 
documentation  of  the  ambient  conditions  for  the  experiment,  or  to  modeling  deficiencies.  Our  aim  in  the 
present  section  is  to  highlight  some  of  the  reasons  for  the  uncertainties  surrounding  deposition  values. 

Gaseous  and  particulate  removal  rates  from  a  polluted  atmosphe»-'  *o  environmental  surfaces 
historically  have  been  reported  as  deposition  velocities,  vj,  defined  ;  ie  flux  to  the  surface  wc|0 
divided  by  the  airborne  concentration,  C,  i.e.. 


vd  *  -wc|0/C(z) 


(6.4.1) 


Since  a  concentration  gradient  is  necessary  to  deliver  a  flux  to  the  surface  (except  in  the  nominal  case 
of  simple  gravitational  settling  of  large  particles),  vj  is  necessarily  a  function  of  the  height  at  which 
the  normalizing  concentration  is  measured.  Figure  6.4.1  gives  the  deposition  velocities  tabulated  by 
Sehmel  (1980)  for  a  number  of  different  surfaces  as  measured  for  one  particular  species.  Dry  deposition 
involves  such  diverse  phenomena  as  leaf  chemistry,  plant  structure,  land  use  as  well  as  all  the  factors 
involved  in  the  turbulent  transport  of  mass  and  momentum  discussed  in  Chapters  1  through  5.  It  is  more 
convenient  to  deal  with  the  resistance  of  the  surface  layer  to  deposition  which  is  defined  as  the 
reciprocal  of  the  deposition  velocity, 


Rd  * 
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This  concept  has  the  advantage  that  it  can  conveniently  be  divided  into  the  resistance  of  the  various 
layers  through  which  the  species  of  interest  must  pass  to  reach  the  surface. 
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Figure  6.4.1  -  Iodine  deposition  velocity  summary  as  given  by  Sehmel  (1980).  Reference 
nunber  refers  to  Sehmel' s  tables. 


When  the  boundary  layer  Is  divided  into  four  different  regions  as  shown  In  Figure  6.4.2,  the  total 
resistance  to  deposition  Is  the  sum  of  resistance  presented  by  each  separate  region  plus  a  fifth 
resistance  due  to  absorption  or  reaction  at  the  surface. 


Vd  ■  (Rt  ♦  ♦  R,  ♦  R4  ♦  Rs)*1  (6,4.3) 
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Figure  6.4.2  -  Four  different  regions  of  the 
planetary  boundary  layer  for  study  of  the  aerodynamic 
resistance  to  SO^  deposition. 


For  flow  over  a  i.  surface.  Region  3,  the  canopy,  does  not  exist  and  R3  •  0.  When  a  vegetal  canopy 
exists  at  the  s  ~e,  the  effective  surface  area  is  increased,  and  R„  Is  decreased  and  may  be 
incorporated  within  Ra.  A  fifth  resistance,  a  surface  chemical  or  biological  resistance,  must  be  added  to 
the  four  aerodynamic  resistances  in  the  most  general  case.  Each  of  these  layers  can  be  considered 
separately. 

The  contribution,  Ru  from  Region  1,  the  outer  boundary  layer,  can  generally  be  neglected  because  the 
height  at  which  the  deposition  velocity  is  generally  desired  falls  within  the  surface  layer  covered  by 
Region  2.  If  the  contribution  is  desired,  then  the  methods  discussed  in  Section  6.1  need  to  be  applied. 

The  resistance  for  Region  2  can  be  obtained  directly  from  discussed  in  Section  3.1,  since  the 
normalized  species  gradient  is  assumed  to  be  the  same  function  of  z/L  as  the  normalized  temperature 
gradient.  This  assumption  is  based  on  the  similarity  of  the  fundamental  heat  and  species  diffusion 
equations.  From  Eqs.  (6.4.1  and  6.4.2) 
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This  Integral  was  given  in  Eq.  (3.1.4)  for  unstable  conditions  and  in  Eq.  (3.1.11)  for  stable  conditions. 
For  neutral  conditions  it  reduces  simply  to 
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The  departure  of  the  resistance  as  a  function  of  height  and  stability  from  the  log  variation  is  given  in 
Figure  6.4.3.  The  figure  shows  a  strong  effect  of  both  stability  and  height  on  deposition  during  stable 
conditions  and  very  little  effect  during  unstable  conditions. 

The  smooth  surface  resistance  is  determined  by  the  Schmidt  number,  v/K c,  of  the  species  diffusing  in 
air,  varying  approximately  as  ( v/Kr)*/3.  For  real  surfaces  covered  with  some  form  of  vegetation,  this 
must  be  combined  with  a  canopy  model,  such  as  discussed  in  Section  3.6,  to  model  the  resistance  to 
deposition  within  the  canopy.  Figure  6,4.4  shows  the  sensitivity  of  the  resistance  of  a  typical  forest 
canopy  with  respect  to  three  Important  variables  according  to  the  model  presented  by  tewellen  and  Sheng 
(1980a).  Everything  else  being  equal,  the  canopy  resistance  to  deposition  decreases  as  the  laminar 
dlffusivlty  of  the  depositing  species  Increases,  as  the  ratio  of  the  plant  leaves'  wetted  area  to  the 
frontal  area  projected  into  the  wind  Increases,  and  also  with  increasing  wind  speed. 
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Figure  6.4.3  -  Incremental  resistance  to  deposition 
as  a  function  of  height  and  stability  as  given  by 
surface  layer  functions. 


y.(em/HC 


0/v  "Sc"' 

Figure  6.4.4  -  Gaseous  deposition  velocity  at  the  top 
of  the  summer  forest  as  a  function  of  three  variables 
(lewellen  and  Sheng,  1980a).  Plant  area  index  *  3. 
wet  foliage,  neutral  stability. 

(a)  As  a  function  of  Schmidt  number  (u^h  “  5  m/sec, 
Ay,/Af  *  10,  u*  »  70  cm/sec). 
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Aw  /A) 


Figure  6.4.4  -  (b)  As  a  function  of  leaf  wetted  area 
to  frontal  area  (D/v  =  1,  U2h  =  6  m/sec). 


Figure  6.4.4  -  (c)  As  a  function  of  U2h  (0/ v  «  l, 
VAf  *  10). 


Different  component  resistances  will  be  Important  In  Eq.  (6.4.3)  under  different  conditions. 
Representative  estimates  for  some  of  the  terms  under  various  conditions  are  given  in  Table  6.4.1.  There 
are  certain  conditions  where  any  -ne  of  the  Individual  resistances  may  be  the  dominant  one.  In  view  of 
Eq.  (6.4.3)  and  Table  6.4.1,  it  should  not  be  surprising  that  such  a  large  scatter  In  numerical  values  of 
deposition  velocity  exists  as  given  In  Figure  6.4.1  for  a  single  species.  It  Is  also  clear  that  models 
which  parameterize  deposition  velocity  as  constant  for  a  particular  species  are  subject  to  large  errors. 
An  attempt  at  a  more  detailed  parameterization  Is  given  by  Lewellen  and  Sheng  (1980a). 
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Table  6.4.1.  Representative  values  of  deposition  resistance 
terms  In  Eq.  4.4.3. 

(sec/cm) 


Minimum 

Low 

Typical 

HI  oh 

Maximum 

(1)  Aerodynamic 
turbulent 

0.05 

0.2 

1 

r 

100 

(R,  +R2) 

Unstable 

Modest 

Light 

Stable  winds 

Light 

flow  over 

neutral 

neutral 

stable 

rough 

surface 

winds 

winds 

winds 

(2)  Laminar 

sublayer 

0.01 

0.1 

1 

20 

1000 

(r3  ♦  r4) 

Gas  species 

Gas  species 

0.01  urn 

0.5  wm 

0.5  urn 

transferred 

transferred 

particles 

particles 

particles 

to  pine 

to  smooth 

on  pine 

deposited  on 

deposited 

forest 

surface 

forest 

pine  forest 

on  smooth 

under  moderate 
wind  conditions 

surface 

(3)  Surface  Rc 
(for  S02) 

0.02 

0.2 

0.5 

1 

10 

(Water) 

(Calcarious 

(Grass  or 

(Pine  with 

(Pine 

soil) 

wheat) 

open  stomata) 

with 

closed 

stomata) 
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VII.  Concluding  Remarks 


Rapid  improvements  in  techniques  for  predicting  the  dominant  processes  in  the  lower  atmosphere  may  be 
expected  over  the  next  few  years  due  to  the  fortunate  combination  of  the  availability  of  improved 
instrunentation  for  probing  this  layer,  the  availability  of  turbulent  transport  models  for  analyzing  the 
sensitivity  of  turbulent  transport  in  this  region  to  particular  variables,  and  the  probable  availability 
of  continued  funding  due  to  the  recognition  of  the  importance  of  transport  in  this  layer  to  a  number  of 
engineering  applications. 

A  variety  of  remote  sensors  based  on  radio,  optical  or  acoustic  techniques  are  currently  being 
tested.  As  some  of  these  prove  their  ability  to  give  good  spatial  and  temporal  resolution  at  a  relatively 
inexpensive  unit  cost  per  observation,  new  applications  will  be  found.  The  Phoenix  Project  (Hooke,  1980) 
is  an  example  of  the  type  of  projects  which  may  be  expected  to  play  a  key  role  in  this  development.  The 
interaction  of  modeler  and  observer  is  important  to  both.  Although  it  is  generally  recognized  that  the 
modeler  must  have  empirical  inputs  for  model  development  and  evaluation,  it  is  equally  true  that  the  giant 
data  bases  generated  by  observational  programs  need  inputs  from  models  for  guidance  in  determining  which 
parameters  are  most  critically  needed.  It  is  only  slightly  wastefull  when  extraneous  information  is 
collected  and  filed,  but  much  more  inefficient  when  the  absence  of  an  important  variable  severely  limits 
the  utility  of  a  large  data  set.  Progress  is  most  efficient  when  there  is  continued  iteration  between 
model  formulation  and  observational  testing  to  determine  the  model's  deficiencies  and  supply  the  empirical 
information  needed  for  improved  model  parameterization. 

Some  of  the  major  applications  expected  to  provide  the  support  to  drive  the  continued  development  of 
transport  modeling  in  the  lower  atmosphere  have  been  discussed  in  some  detail  in  the  previous  chapter. 
There  are  many  other  areas  which  may  contribute  to  this  support,  since  they  will  certainly  benefit  from 
the  expanded  capabilities.  These  run  the  spectrum  from  agricultural  planning  to  zeppelin  operations. 
Some  areas  which  may  be  expected  to  have  a  particularly  strong  interaction  with  our  increasing  knowledge 
about  turbulent  transport  in  the  lowest  1  km  of  the  atmosphere  are  represented  in  the  two  volume 
Proceedings  of  the  5th  International  Conference  on  Wind  Engineering  (Cermak,  1980).  Another  area  where 
improved  understanding  of  the  boundary  layer  should  reap  substantial  dividends  is  in  the  parameterization 
of  surface  influences  in  large-scale  weather  forecasting  models  (Driedonks  and  Tennekes,  1981).  Perhaps 
the  strongest  impetus  to  improved  prediction  techniques  may  come  from  the  military.  Not  only  is  there  the 
direct  requirement  of  forecasting  the  natural  micrometeorological  environment  in  which  men  and  weapons  are 
called  to  operate,  but  there  is  the  potential  use  of  tactical  modifications  to  this  environment. 

Many  important  problems  fall  under  the  area  of  the  interaction  between  energy  conversion  and  the 
environment.  The  dispersal  of  the  residue  from  a  coal-fired  power  plant  is  an  obvious  example.  Less 
obvious  are  such  things  as  the  impact  of  the  accuracy  of  dispersal  models  on  emergency  planning  around 
nuclear  power  plants;  the  possible  adverse  impact  of  cooling  towers,  or  ponds,  on  local  fog  episodes; 
and  the  impact  of  improved  estimates  of  surface  transfer  coefficients  on  questions  surrounding  the  global 
atmospheric  buildup  of  CO^  as  a  limit  to  our  use  of  fossil  fuels.  The  strong  interaction  extends  into  the 
areas  of  renewable  energy  sources  and  energy  conservation.  The  choice  of  windmill  sites  and  modes  of 
operation  depends  on  knowledge  about  the  wind  behavior  close  to  the  surface.  The  maintenance  of  active 
solar  energy  conversion  systems  depends  on  wind  damage  and  dust  erosion.  Architects  are  recognizing  the 
influence  of  siting  and  landscaping  on  heating  and  cooling  requirements  for  buildings. 

Many  applications  could  profit  from  a  regional  boundary- layer  meteorological  model  based  on  inputs 
from  regular  weather  forecasts  along  with  local  topographic  features,  and  local  sources  and  sinks  of 
thermal  energy,  hunidity,  and  effluents  of  particular  interest.  Some  research  models  currently  are 
attempting  this  (e.g.,  Pielke,  1974;  Tapp  and  White,  1976;  and  Anthes  and  Warner,  1978).  One  model 
(Long  and  Shaffer,  1975)  attempts  this  on  an  operational  basis,  but  on  a  relatively  large  grid  (80  km 
between  horizontal  grid  points)  so  that  micrometeorological  features  cannot  be  resolved.  Advancements  in 
three  specific  areas  are  needed  to  make  such  a  regional  model  an  attractive  possibility.  First,  efficient 
algorithms  for  providing  sufficient  physical  constraints  must  be  developed  to  allow  relatively  fine  mesh 
calculations  at  reasonable  costs.  Second,  and  perhaps  most  difficult,  is  the  problem  of  supplying 
boundary  conditions  with  adequate  resolution  in  time  and  space  to  drive  the  model.  Finally,  there  is  the 
requirement  of  processing  the  large  flow  of  output  data  in  a  useful  manner  for  prospective  users. 
Continued  reductions  In  unit  computing  costs  contribute  to  a  high  probability  of  solving  the  first  and 
third  of  these,  but  there  is  no  assurance  as  to  how  well  the  second  can  be  solved.  The  challenge  of 
following  the  flow  In  real  time  is  to  provide  select  data  such  as  might  be  available  from  meteorological 
stations  within  the  domain  or  satellite  remote  sensing  to  supplement  the  sparse  boundary  data  as  proper 
constraints  on  the  flow.  The  forecast  capability  Involves  this  plus  the  basic  question  of  the  limits  to 
predictability  in  a  turbulent  atmosphere  as  discussed  by  Haltiner  and  Williams  (1980). 
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14.  Abstract 

The  interdependence  between  the  turbulent  transport  of  mass,  momentum,  and  energy 
through  the  lower  levels  of  the  atmosphere,  and  the  distribution  of  wind,  temperature,  and 
species  within  this  layer  are  reviewed.  Particular  emphasis  is  placed  on  models,  either 
analytical  or  numerical,  which  have  a  basic  theoretical  foundation  in  turbulent  transport 
modeling.  The  results  of  example  model  solutions  are  used  to  discuss  such  micro- 
meteorological  problems  as:  wind  shear  and  turbulence  around  airports;  the  prediction  of 
low-level  clouds  and  fog;  the  dispersion  of  industrial  pollutants;  the  interaction  of  electro¬ 
magnetic  radiation  with  turbulent  fluctuations  in  humidity  and  temperature ;  and  the 
evolution  of  organized  features  within  the  boundary  layer. 
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